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FOREWORD 


This  report  documents  computer  software  that  was  developed  as  part  of  an 
investigation  of  the  non-isothermal  viscoplastic  response  of  the  titanium  alloy 
Timetal®21S.  The  software  was  completed  and  tested  by  Joseph  L.  Kroupa  of  the 
Structural  Integrity  Division,  University  of  Dayton  Research  Institute,  Dayton, 
OH.  The  research  was  conducted  at  the  Materials  Behavior  Branch,  Metals  and 
Ceramics  Division,  Materials  Directorate,  Wright  Laboratory  (WL/MLLN) 
Wright-Patterson  AFB,  OH,  tmder  Contract  No.  F33615-94-C-5200.  The  contract 
was  administered  tmder  the  direction  of  WL/MLLN  by  Mr.  Jay  R.  Jira.  Dr.  Noel 
E.  Ashbaugh  was  the  Principal  Investigator  and  Dr.  Joseph  P.  Gallagher  was  the 
Program  Manager. 


THIS  SOFTWARE  AND  ANY  ACCOMPANYING  DOCUMENTATION 
IS  RELEASED  “AS  IS.”  THE  U.S.  GOVERNMENT  MAKES  NO 
WARRANTY  OF  ANY  KIND,  EXPRESS  OR  IMPLIED,  CONCERNING 
THIS  SOFTWARE  AND  ANY  ACCOMPANYING  DOCUMENTATION, 
INCLUDING,  WITHOUT  LIMITATION,  ANY  WARRANTIES  OF 
MERCHANTABILITY  OR  FITNESS  FOR  A  PARTICULAR  PURPOSE.  IN 
NO  EVENT  WILL  THE  U.S.  GOVERNMENT  BE  LIABLE  FOR  ANY 
DAMAGES,  INCLUDING  ANY  LOST  PROFITS,  LOST  SAVINGS  OR 
OTHER  INCIDENTAL  OR  CONSEQUENTIAL  DAMAGES  ARISING  OUT 
OF  THE  USE,  OR  INABILITY  TO  USE  THIS  SOFTWARE  OR  ANY 
ACCOMPANYING  DOCUMENTATION,  EVEN  IF  INFORMED  IN 
ADVANCE  OF  THE  POSSIBILITY  OF  SUCH  DAMAGES. 


vi 


SECTION  1.0 


INTRODUCTION 

Many  research  and  design  projects  involving  structural  response  at 
elevated  temperatures  require  a  visco-plastic  model  to  adequately  describe 
strain-rate  sensitivity  and  time-dependent  behavior  of  material.  In  the  case  of 
Timetal®21S,  a  unified-inelastic-strain  theory  consisting  of  a  modified  flow  rule 
of  a  Bodner-Partom  [1]*  works  well  to  describe  its  visco-plastic  response.  In 
addition  to  sophisticated  material  models,  many  design  projects  also  require 
stress-strain  analysis  in  two  and  three-dimensions  (e.g.,  multiaxial  stress  states) 
with  complex  boundary  conditions.  This  is  easily  obtained  with  the  finite 
element  method.  Thus,  the  integration  of  the  Bodner-Partom  model  into  a  well 
established  finite  element  package  makes  for  a  very  powerful  tool  in  solving 
problems  that  require  visco-plastic  material  response  and  multiaxial  stress  states. 
Due  to  the  advanced  nature  of  these  visco-plastic  theories,  many  finite  element 
programs  do  not  include  these  formulations  directly  into  their  codes.  However, 
most  commercially  available  finite  element  programs,  such  as  ABAQUS  [2],  do 
allow  the  incorporation  of  user-defined  subroutines  that  would  integrate  the 
visco-plastic  theories  into  their  algorithms.  This  manual  describes  a  imified 
inelastic  strain  theory  for  Timetal®21S  and  provides  the  associated  user-defined 
subroutines  that  integrate  the  theory  into  ABAQUS.  The  manual  is  designed  to 
assist  the  designer  and  researcher  in  using  and  understanding  these  subroutines 
within  ABAQUS. 
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This  document  shares  several  similarities  with  a  predecessor  manual  that 
describes  the  implementation  and  use  of  the  Bodner-Partom  model  in  ADINA 
[3].  One  similarity  is  the  description  of  the  visco-plastic  theory  which  appears 
both  manuals  for  completeness.  The  second  similarity  is  that  the  same  numerical 
algorithms  are  used.  The  most  obvious  difference  is  that  the  present  document 
describes  the  interface  details  associated  with  ABAQUS.  In  addition,  this  manual 
includes  more  information  on  error  debugging  and  common  usage  errors. 

In  conjimction  with  improvements  to  the  user  manual,  improvements  to 
the  numerical  algorithms  have  also  been  completed.  One  major  improvement  is 
the  coding  of  a  more  generalized  form  of  the  Bodner-Partom  equations,  which 
allows  more  temperature-dependent  material  properties  than  previously 
considered.  In  association  with  this  modification,  a  second  set  of  material 
constants  is  provided  that  better  capture  the  creep  behavior  of  Tunetal®21S  than 
those  proposed  by  Neu  [4].  The  resultant  set  of  algorithms  can  readily  handle 
any  temperature  dependent  material  parameters  for  a  wide  variety  of  materials. 

Other  improvements  include:  1)  a  unified  set  of  subroutines  for  both  the 
two-  and  three-dimensional  elements;  2)  new  algorithms  that  accommodate  the 
plane  stress  case  without  further  modifications;  3)  a  more  straight  forward 
iterative  loop;  and  4)  a  special  flag  to  automatically  reduce  the  time-step  size. 
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SECTION  2.0 


BODNER-PARTOM  CONSTITUTIVE  MODEL 


This  section  presents  the  Bodner-Partom  constitutive  equations  and  the 
material  parameters  that  describe  the  Timetal®21S  response.  Several 
formulations  of  the  Bodner-Partom  model  can  characterize  inelastic  deformation 
under  a  variety  of  conditions,  including  anisotropic,  isothermal,  and/or  non- 
isothermal  material  response.  The  formulation  presented  in  the  Section  2.1 
parallels  the  isotropic  nonisothermal  theory  of  Chan,  Bodner,  and  Lindholm 
[5,6,7].  The  terminology  is  similar  to  that  of  Chan  and  Lindholm  [8].  Section  2.2 
contains  model  parameters  that  characterize  the  Timetal®21S  material. 

2.1  Theory 


The  first  assumption  in  this  analysis  is  the  decomposition  of  the  total 
strain,  eJ®*  ,  into  elastic,  thermal,  and  inelastic  components.  This  decomposition 
is  expressed  as: 


^  el  ,  th  .  in 


Eq.  1 


The  elastic  strains,  ,  depend  on  the  current  stress  state,  elastic  modulus  E,  and 
Poisson's  ratio  v.  The  thermal  strain  components,  Ey'’,  equal  the  product  of  the 
coefficient  of  thermal  expansion  and  the  difference  between  the  current  and 
reference  temperatures.  The  Bodner-Partom  flow  rule  governs  the  evolution  of 
the  inelastic  strains,  e|"  . 
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As  opposed  to  other  unified  inelastic  strain  formulations,  this  theory 
describes  the  directional  (kinematic)  hardening  with  a  special  directional 
hardening  term.  Other  theories  represent  directional  hardening  phenomena 
with  a  "back-stress"  modified  equivalent  stress  [9,10,11].  Introduction  of  the 
directional  hardening  term  alters  the  Bodner-Partom  flow  rule  [1]  by  replacing 
the  previously  known  variable  "drag-stress"  with  the  sum  of  isotropic  and 
directional  hardening  terms,  Z'  and  7?,  respectively.  These  two  hardening  terms 
enter  into  the  inelastic  strain  rate  equation  or  flow  law: 


if  iz' 


Eq.2 


where  Dj,  is  the  limiting  strain  rate,  Sj^are  the  components  of  deviatoric  stress,  and 


The  evolution  of  Z’  and  Z°  have  similar  empirical  forms.  Each  equation 
consists  of  a  hardening  term,  a  thermal  recovery  term,  and  a  temperature  rate 
term.  The  isotropic  hardening  evolution  equation  with  these  three  terms  is: 

fz'-Z,Y’  (z.-z'^sz,  . 

Z'=m,(z,-Z')W"-A,Z,  Eq.3 


Eq.3 


where  the  inelastic  work  rate  is  given  by: 


<7.e'" 
y  y 


Eq.4 


The  initial  value  of  the  isotropic  hardening,  Z\0),  is  Z^.  The  material  parameters 

associated  with  the  isotropic  hardening  evolution  are  mj,  Zg,  Z^,  Z^,  Aj,  and  r^. 

dZ, 

The  thermal  differential  term  appropriately  scales  the  isotropic  hardening 
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variable  when  inelastic  deformation  and  thermal  recovery  occur  imder 
nonisothermal  conditions;  th\is’  Z’  is  prevented  from  passing  through  maximum 
or  minimum  values  with  temperature  changes.  The  treatment  of  this  thermal 
differential  term  is  more  consistent  with  the  work  of  McDowell  [12]  and  others 
[13,14,15,16]  than  those  proposed  by  Chan,  Bodner  and  Lindholm  [5]. 

The  scalar  product  of  a  state  variable,  and  a  imit  stress  vector,  Uj^, 
yields  the  magnitude  of  the  directional  hardening  term: 

Eq,5 

where: 


U..  —  •—===■ 

and: 

_ ^ 

^/AA■ 

The  initial  directional  hardening  variable,  Z^{0)  is  set  to  zero.  The  material 

constants  associated  with  the  directional  hardening  evolution  equation  are  m^, 

dZ, 

Z3,  A^,  and  rj.  The  temperature  differential  term  appropriately  scales  the 

directional  hardening  variable  when  inelastic  deformation  and  thermal  recovery 
occurs  imder  nonisothermal  conditions.  In  particular,  without  these  differential 


Eq.7 


Eq.8 


i3;j  =  m2(Z3Uij-/3ij)W“-A2Zj, 


V^kl^kl 


Pi)  aZ3 


t,  Eq.6 
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terms  the  directional  hardening  accrued  at  one  temperature  may  exceed  the 
limiting  value  Zj  at  another  temperature,  which  is  not  physically  realistic. 

Table  2-A  summarizes  the  15  material  parameters  that  characterize  the 
strain-rate  sensitivity  and  time-dependent  behavior  of  Timetal®21S  [17].  The 
number  of  material  parameters  is  effectively  less  than  those  in  Table  2-A  after 
applying  the  usual  assumptions;  Aj  =  A^,  r,  =  r^,  and  =  Z^.  The  number  of 
temperature  dependent  constants  is  minimized  to  reduce  the  amoimt  of  required 
data.  Other  materials,  such  as  revised  Timetal®21S  constants  [18],  Mar-M47  [5] 
or  B1900+Hf  [6],  require  different  temperature-dependent  parameters.  The 
numerical  algorithms  provided  in  Appendix  1  assume  that  all  the  material 
properties  are  temperature  dependent,  except  and  Zj.  The  algorithms  also 
assume  that  all  material  properties  are  unique,  i.e.,  Aj  does  not  necessarily  equal 


Table  2-A 


Material  Parameters  in  Bodner-Partom  Constitutive  Model 


Parameter 

Units 

Description 

E 

MPa 

Elastic  modulus 

V 

— 

Poisson’s  ratio 

a 

(°cr 

Coefficient  of  tiiermal  expansion 

D.. 

s' 

Limiting  shear  strain  rate 

z. 

MPa 

Initial  value  of  the  isotropic  hardening  variable 

z, 

MPa 

Limiting  (maximum)  value  of  Z 

z. 

MPa 

Fully  recovered  (minimum)  value  of  Z' 

Z3 

MPa 

Limiting  (maximum)  value  of  Z° 

m, 

(Mpa)' 

Hardening  rate  coefficient  of  Z' 

m. 

(Mpa)' 

Hardening  rate  coefficient  of 

n 

— 

Strain  rate  sensitivity  parameter 

A, 

s" 

Recovery  coefficient  for  Z‘ 

A, 

s’ 

Recovery  coefficient  for  Z^ 

r, 

— 

Recovery  exponent  for  Z' 

r. 

— 

Recovery  exponent  for 
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2.2  Timetal®21S  Material  Parameters 


The  material  parameters  for  the  Bodner-Partom  model  with  directional 
hardening  were  determined  from  Timetal®21S  using  monotonic,  cyclic,  and 
creep  test  data  [17].  The  parameters  are  valid  for  a  wide  range  of  strain  rates  (10'^ 
to  10'^  s'^)  and  temperatures  (23°  -  815°C).  The  strategy  for  determming  the 
parameters  involves  a  number  of  steps.  First,  the  temperahire-dependent 
parameters  are  identified  and  the  values  of  the  temperature-independent 
parameters  are  estimated.  Then  at  each  temperature  where  experimental  data 
are  available,  the  temperature-dependent  parameters  are  determined  through  an 
iterative  process  to  minimize  the  differences  between  model  simulations  and 
experiments.  Mathematica  [19]  is  used  to  generate  the  model  simulations. 
Similar  to  other  inelastic  strain  theories,  the  set  of  material  parameter  for  any 
particular  load  case  is  not  imique.  Thus,  for  a  given  set  of  experimental  load 
responses,  a  range  of  values  is  suitable  for  each  material  parameter. 

The  resulting  set  of  temperature-dependent  parameters  becomes 
continuous  with  temperature  as  the  range  of  possible  values  is  decreased  for 
each  parameter.  The  response  can  be  very  sensitive  to  small  changes  in  certain 
material  parameters  with  temperature,  especially  in  the  transition  regimes 
between  different  deformation  mechanisms.  For  Timetal®21S,  a  transition 
region  for  inelastic  behavior  occurs  between  482°C  for  plasticity  and  at  650°C  for 
power  law  creep.  Thus,  an  anomalous  change  in  the  saturated  stress  level  can 
occur  if  linear  interpolation  of  material  parameters  is  used  within  this  transition 
region. 
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Since  no  experimental  data  are  available  within  this  transition  region, 
values  for  temperature-dependent  parameters  are  chosen  between  482°C  and 
650°C.  Thus,  the  resulting  saturated  stress  is  smooth  and  continuous  with 
temperature,  reducing  the  anomalous  effects  of  linear  interpolation.  The  final 
version  of  the  material  parameters  for  Timetal(E)21S  is  shown  in  Table  2-B  [4]. 

Note  that  at  several  temperatures  Z3  is  larger  than  Z^.  This  introduces  a 
potential  convergence  problem  upon  load  reversals  at  high  strain  rates.  A 
revised  series  of  material  parameters  for  Timetal®21S  [18],  provided  in  Table  2- 
C,  is  proposed  to  eliminate  this  problem,  tn  addition,  the  revised  parameters  are 
calibrated  to  better  capture  lower-stress  creep  response  than  the  previous 
constants,  which  are  fit  to  the  widest  range  of  stress  levels  and  temperatures 
possible. 
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Table  2-B 


Bodner-Partom  Material  Parameters  for  Timetal®21S  141 

Temperature-Independent  Constants 

m,  =  0  MPa  '  (no  isotropic  hardening)  Z,  =  1600  MPa 
r,  =  r^  =  3  D„  =  1  X 10^  s’ 

V  =  0.34 


Temperature-Dependent  Constants 


Temp. 

°C 

E 

MPa 

a* 

10^ 

(°C)-’ 

n 

z.  =  z, 

MPa 

Z3 

MPa 

m, 

(Mpa)' 

> 

II 

23 

112000 

6.31 

4.8 

1550 

100 

0.35 

0 

260 

108000 

7.26 

3.5 

1300 

300 

0.35 

0 

315 

<t' 

<t> 

0 

0 

390 

0 

<1> 

365 

«t> 

<1> 

0 

0 

500 

0 

<1> 

415 

<(> 

0 

0 

660 

0 

465 

<!> 

<!> 

0 

0 

960 

0 

<1> 

482 

98100 

8.15 

1.7 

1100 

1100 

5 

0.0076 

•500 

<}> 

2.5 

0 

1300 

0 

<t> 

525 

<l> 

<!> 

1.28 

0 

1670 

0 

<l> 

550 

<|) 

1.1 

0 

2100 

0 

<1> 

575 

<|) 

0.97 

0 

2600 

0 

<t> 

600 

(t* 

0.82 

0 

3700 

10 

650 

86600 

8.83 

0.74 

1000 

3800 

10 

0.21 

760 

77200 

9.27 

0.58 

600 

4000 

15 

1.0 

815 

72000 

9.49 

0.55 

300 

4100 

30 

2.0 

*  Secant  a  with  T„  =  25°C 

0  -  Linear  interpolate  between  values  given  in  table. 

(|)  -  Use  functions  to  determine  values. 

Bold  were  values  determined  based  on  experiments. 

Italics  are  values  that  describe  smooth  and  continuous  saturated  stress  change. 
Functions 


A,  =  C5,  exp 
A,  =  C5,  exp 


f 

V 

r 


'61 


T  +  273 


^62 

T  +  273 


,  with  C5,  = 
,  with  €52  = 


5.8  X 10^  s' ,  C„ 
5.8  X 10^  s  ’ ,  C„ 


1.37xlO^CandTin°C 

1.37xlO‘C,andTin°C 


10 


Table  2-C 


Revised  Bodner-Partom  Material  Parameters  for  Timetal®21S  1181 


Temperature-Independent  Constants 

Z,  =  3200MPa 
D„=lxl0^s-' 
Z3  =  180  MPa 


Temperature-Dependent  Constants 


‘•r — = — =- 

Temp. 

°C 

E 

MPa 

a* 

10-*  (°C)' 

n 

=  z. 

MPa 

m, 

(Mpa)' 

A,  =  . 
s-' 

23 

112000 

6.31 

1.58 

2950 

0.5 

0 

260 

108000 

7.26 

1.35 

2650 

0.5 

0 

566 

92700 

8.49 

1.09 

800 

2.5 

0.3 

650 

86600 

8.83 

0.94 

200 

7.0 

2.0 

760 

77200 

9.27 

0.575 

190 

37.0 

4.5 

815 

72000 

9.49 

0.55 

150 

50 

7.5 

lUj  =  4  MPa' 
r,  =  Tj  =  3.5 
V  =  0.34 


*Secant  a  with  T„  =  25'’C 


SECTION  3.0 


SUBROUTINE  OPERATION 

3.1  Subroutine  Installation  and  Execution 

The  installation  of  the  user-defined  subroutines  into  ABAQUS  is  fairly 
straight  forward.  The  user-subroutines  are  compiled  and  linked  automatically 
when  executing  the  code.  Thus,  a  permanent  version  of  executable  code  is  never 
created. 

The  user-subroutines  can  be  incorporated  into  ABAQUS  in  two  manners. 
One  method  is  to  include  the  subroutines  directly  into  the  ABAQUS  input  file 
after  a  *USERSUBROUTlNE  card.  This  method  works  well  for  one  time  use,  but 
is  not  best  suited  for  multiple  input  files.  For  use  with  several  input  files,  the  best 
method  is  to  include  the  file  name  in  the  ABAQUS  command  statement.  For 
example,  if  the  user-subroutines  reside  in  the  file  usub_db .  f ,  then  the  command 
statement  for  ABAQUS  would  be: 

abaqus  job=jobnaine  user=usub_db 

Execution  of  this  command  will  look  for  the  input  files  j  obname .  inp  and 
usub_db .  f ,  This  command  will  automatically  compile  and  link  the  user- 
subroutines  in  usub_db .  f  to  create  an  executable  version  of  ABAQUS.  The 
program  then  proceeds  to  solve  the  problem  described  in  j  obname .  inp . 
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3.2  Material  Property  Input 


Different  from  ADINA,  where  temperature  dependent  properties  are 
assumed,  ABAQUS  assumes  nothing  special  about  the  material  properties.  As  a 
result,  the  properties  are  read  in  as  a  single  array,  which  requires  a  special 
subroutine  to  extract  the  data.  Figure  3.1  gives  an  example  how  the  Bodner- 
Partom  material  constants  of  Table  2-B  appear  in  the  ABAQUS  input  file. 


♦MATERIAL,  NAME=B21S_3 
♦EXPANSION,  ZERO=900.0 
9.7787E-6  23.0 

1.0713E-5  260.0 


More  cards  here 

1.2323E-5  760.0 
1.2467E-5  815.0 
1.2689E-5  •  900.0 


♦USER  MATERIAL,  C0NSTANTS=176 


0 

13. 

16. 

1.0E4 

1600.0 

5.8E5 

1.37E4 

5.8E5 

1.37E4 

2.0 

10.0 

0.75 

l.E-3 

0.0 

0.0 

0.0 

23.0 

260.0 

315.0 

365.0 

415.0 

465.0 

482.0 

500.0 

525.0 

550.0 

575.0 

600.0 

650.0 

760.0 

815.0 

900.0 

•112000.0 

108030.0 

106130.0 

104090.0 

101740.0 

99085.0 

98113.0 

97045.0 

95497.0 

93873.0 

92172.0 

90395.0 

86612.0 

77000.0 

72000.0 

63000.0 

0.340 

0.340 

0.340 

0.340 

0.340 

0.340 

0.340 

0.340 

0.340 

0.340 

0.340 

0.340 

0.340 

0.340 

0.340 

0.340 

4.80 

3.50 

3.054 

2.649 

2.243 

1.838 

1.700 

1.500 

1.280 

1.100 

0.970 

0.820 

0.740 

0.580 

0.550 

0.550 

1550.0 

1300.0 

1250.4 

1205.4 

1160.4 

1115.3 

1100.0 

1089.3 

1074.4 

1059.5 

1044.6 

1029.8 

1000.0 

600.0 

300.0 

300.0 

100.0 

300.0 

390.0 

500.0 

660.0 

960.0 

1100.0 

1300.0 

1670.0 

2100.0 

2600.0 

3700.0 

3800.0 

4000.0 

4100.0 

4300.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.350 

0.350 

1.502 

2.549 

3.597 

4.644 

5.000 

5.763 

6.822 

7.881 

8.941 

10.00 

10.0 

15.0 

30.00 

30.00 

3.0 

3.0 

3.0 

3.0 

3.0 

3.0 

3.0 

3.0 

3.0 

3.0 

3.0 

3.0 

3.0 

3.0 

3.0 

3.0 

3.0 

3.0 

3.0 

3.0 

3.0 

3.0 

3.0 

3.0 

3.0 

3.0 

3.0 

3.0 

3.0 

3.0 

3.0 

3.0 

♦DEPVAR 

20 

♦INITIAL  CCNDinaNS,  TYPE=SOLUITCN,  USER 


Figure  3.1  Listing  of  Material  Properties  within  ABAQUS  input  file. 
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From  Figure  3.1,  note  that  the  thermal  expansion  data  are  input  in  the 
same  manner  as  any  other  model  used  in  ABAQUS.  On  the  *USER  MATERIAL 
card,  CONSTANTS=204  indicates  that  204  constants  are  read  in  from  the 
remaining  cards  in  this  data  block  which  are  described  below.  The  *DEPVAR 
card  and  the  following  card  indicate  that  20  state  variables  are  allocated  in  the 
analysis.  The  next  card,  *  INITIAL  CONDITIONS,  is  required  to  flag  the  main 
program  to  call  the  subroutine  SDVINI ,  which  initializes  state  variables. 


The  variables  defined  imder  the  *USER  MATERIAL  card  are  best 
described  with  Figure  3.2.  The  first  constant,  I  PLANE  (=1) ,  is  a  flag  to  use  plane 
stress  elements.  The  second  constant  NTIP  is  the  number  of  temperature 
independent  material  properties  and  solution  control  variables.  The  third 
constant  NTSET  is  the  number  of  points  that  define  the  temperature  array. 
Appearing  next  are  temperature-independent  material  properties,  (e.g.,  D^,  Cjj, 
Qj,  ..,)  and  the  solution  control  variables  (e.g.,  XNEWDT,  INTER...).  The 
temperature  array  T{1),T{2),  ...  T (NTSET)  appears  next,  which  is 

preceded  by  the  elastic  modulus  array  TE  { 1 )  ,  TE  ( 2 )  ,  ...  TE  (NTSET) .  The 

following  sets  of  arrays  define  Poisson's  Ratio  (v),  n,  Zj,  Z3,  mj,  m^,  rj,  and  r^ 
respectively. 


TJISERMATEKIALCX:^CrANIS=NIOrAL 


IPLRNE 

NTIP 

NTSET 

DO 

C62 

XNEWDT 

INTER 

RELAX 

T(l) 

T(2) 

T(3) 

.... 

TE(3) 

«... 

.... 

TE  (NTSET) 

....  1NU  (NTSET) 

TN(1) 

IN(2) 

TZ2(1) 

TZ2(2) 

TZ2(3) 

.... 

TZ3(3) 

•  •  .  • 

•  •  *  • 

TZ3  (NTSET) 

....  IMl  (NTSET) 

■ncd) 

TM2(2) 

TRl(l) 

TR1(2) 

TR1(3) 

... 

TR2(3) 

. .  . 

. .  .  1^2  (NTSET) 

Z1 

C51 

C61 

C52 

TOLER 

extra 

extra 

extra 

.... 

T  (NTSET) 

TE{1) 

TE(2) 

TNU(l) 

TNU(2) 

TNU(3) 

.... 

'IN(3) 

.... 

....  TN  (NTSET) 

.  •  .  • 

TZ2  (NTSET) 

TZ3(1) 

TZ3(2) 

TMl(l) 

TM1(2) 

TM1(3) 

«... 

TM2(3) 

.... 

....  TIC  (NTSET) 

.... 

TRl  (NTSET) 

TR2(1) 

TR2(2) 

Figure  3.2  Bodner-Partom  material  property  array  use  with  ABAQUS. 
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3.3  Initializing  State  Variables 


For  the  Bodner-Partom  theory  presented  above,  the  only  state  variable 
that  has  a  non-zero  initial  value  is  Z\  According  to  the  theory,  the  initial  value  of 
7}  is  Zg,  which  is  temperature  dependent.  The  SDVINI  subroutine  within 
ABAQUS,  which  initializes  state  variables,  does  not  automatically  assign  a  value 
of  Z'.  Thus  the  initial  value  of  Z'  must  be  "hard  wired"  into  the  subroutine.  In 
other  words,  the  subroutine  SDVINI  must  be  edited  to  assign  the  appropriate 
values  of  Z'.  The  section  of  code  that  must  be  edited  is  shown  by  the  italic  print 
in  Figure  3.3.  Note  that  the  initial  value  of  Z'  is  assigned  to  state  variable  #13. 
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SUBROUTINE  SDVINI {STATEV,COORDS,NSTATV,NCRDS,NOEL,]SIPr,LAYER,KSPr) 
C 

INCLUDE  'ABA_PARAiyi.INC' 

C 

DIMENSION  STATEV(NSTATV) ,  COORDS  (NCRDS) 
c 

COMMON  /matconst/  E,  G,  a3k,  anu,  dg,  d3k,  dm,  dO, 

&  an,  ami,  am2,  zl,  z2,  z3,  al,  a2,  rl,  r2. 

Sc  dn,  dml,  dm2,  dzl,dz2,dz3,dal,da2,dri,dr2 

c 

statev(l)  =  0.0 
statev(2)  =0.0 
statev(3)  =  0.0 
statev(4)  =0.0 
statev(5)  =0.0 
statev(6)  =0.0 
statev(7)  =0.0 
statev(8)  =0.0 
statevO)  =0.0 
statev(lO)  =0.0 
statev(ll)  =0.0 
statev(12)  =0.0 
c 

c  Initializing  ZI  to  ZO 
c 

c  News  Constants  (Table  2-B)  [4] 
c25 

c25  statev(13)  =  1550. 

c 

c650 

c650  statev(13)  =  1000. 

c 

c900 

c900  statev(13)  =  300. 

c 

c  for  Sander's  Constants  (Table  2-C)  [18] 
c25 

c25  statev(13)  =  2950. 

c650 

statev(13)  =  200. 

c815 

c815  statev(13)  =  150. 

c 

statev{14)  =0.0 
statev(15)  =0.0 
statev(16)  =  0.0 
statev{17)  =0.0 
c 

return 

end 


Figure  3.3  Listing  of  subroutine  SDVINI . 
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The  Bodner-Partom  subroutines  were  extensively  tested  to  minimize 
errors  and  maximize  numerical  efficiency.  Even  though  the  subroutines  were 
thoroughly  tested,  some  conditions  will  result  in  failure  of  the  subroutines  to 
converge  to  a  solution.  Documented  here  are  the  errors  that  were  most 
commonly  encountered  during  this  test  period. 

1.  Wrong  initial  value  of  Z\  This  particular  error  is  most  likely  to  occur  if 
several  analyses  are  required  witib  different  initial  temperatures.  The  user 
should  be  aware  that  if  Z'  is  input  erroneously,  the  numerical  algorithms 
may  fimction  normally  with  erroneous  output.  In  addition,  the 
subroutines  will  not  warn  the  user  if  an  inappropriate  value  has  been 
assigned .  See  Section  3.3  about  assigning  a  value  to  Z' 

2.  Numerical  overflow.  This  error  is  most  likely  to  occur  if  the  time  step  size 
is  too  large  for  the  algorithm.  To  avoid  this  error  assign  INTER  a  larger 
munber.  For  most  cases  INTER=10  is  adequate. 

3.  Nonconvergent  Solutions.  There  is  a  couple  of  reasons  that  would  cause 
this  error  to  occur.  The  most  common  reasons  are  given  here: 

3a.  The  time  step  size  for  the  problem  may  be  too  great.  Increase  the 
value  of  INTER. 

3b.  Discontinuities  in  the  temperature  history  are  likely  to  occur  if 
multiple  time  steps  are  not  used  correctly.  The  default  for  the 
*AMPLIDUTE  CARD  is  for  time  to  start  from  the  beginning  of  the 
solution  step,  rather  than  total  time.  The  flag  "TIME=TOTAL  TIME" 
will  amend  this  problem. 

3c.  A  potential  problem  can  occur  with  stress-reversal  at  fairly  high 
loading  rates  with  the  constants  of  Table  2-B.  This  problem  is  an 
artifact  of  the  chosen  constants  at  650°C  rather  than  the  numerical 
algorithms  (e.g.,  Z3  >  Z^).  Two  solution  to  the  problem  are 
to  limit  the  value  to  Z°  to  be  positive  or  to  use  the  constants  of 
Table  2-C. 


SECTION  4.0 


EXAMPLES  AND  VERIFICATION 

This  section  provides  an  example  of  the  application  of  the  Bodner-Partom 
user-defined  subroutines  with  the  Timetal®21S  material  parameters.  This 
section  presents  several  numerical  exercises  that  verify  the  numerical  accuracy  of 
the  subroutines  in  Appendix  1.  The  differential  equation  solver  within 
Mathematica  [19]  provides  a  very  accurate  base-line  solution  for  the  evaluation  of 
accuracy.  This  section  provides  a  general  survey  of  numerical  accuracy,  while 
Section  5  provides  a  detailed  investigation. 

4.1  Numerical  Examples 

The  file,  easel. inp,  (see  Appendix  2)  is  an  example  ABAQUS  input  file  that 
contains  the  Timetal®21S  material  parameters.  The  example  consists  of  a  single 
axisymmetric  element  imidirectionally  loaded  at  a  constant  strain  rate  of 
833.3xl0'*s'’  and  constant  temperature  of  25°C,  which  has  the  same  conditions  as 
case  1  (discussed  below).  Input  of  the  Timetal®21S  material  properties  starts 
with  the  string  *MATERIAL ,  NAME=B2 1S_3  (see  Appendix  2).  In  these  example 
analyses,  the  thermal  expansion  coefficients  are  set  to  zero  (commented  out) 
since  the  comparisons  only  consider  mechanical  strain  (total  minus  thermal 
strain).  The  output  for  this  numerical  example  is  provided  below,  in  Figure  4.1. 
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4.2  Solution  Verification 


This  section  presents  five  test  cases  that  were  compared  to  the  solution 
obtained  by  Mathematica,  which  has  an  associated  error  of  approximately  10^ 
percent.  The  measure  for  accuracy  is  the  percent  deviation  of  the  finite  element 
stress  solution  from  the  Mathematica  results.  The  percent  stress  deviation,  %D,  is 
written  as: 


%D  =  m%* 


^  ^  ^  math 

^math 


Eq.9 


Where  and  are  the  stresses  obtained  from  the  finite  element  and 

Mathematica  solutions,  respectively.  Smaller  percent  deviation  is  indicative  of  a 
more  accurate  numerical  solution,  rather  than  a  better  fit  with  experimental  data. 

All  cases  simulate  a  monotonically  increasing  tensile  loaded  element  with 
constant  mechanical  strain  rates  of  either  833.3xl0'*s’‘  or  8.33xl0’*s''.  Table  4-A 
summarizes  all  the  test  conditions  and  location  of  the  Mathematica  solution  files 
in  the  appendices.  The  verification  cases  compare  well  within  1.0  %  deviation  for 
both  sets  of  Bodner-Partom  material  constants  from  Tables  2-B  and  2-C.  Results 
are  shown  only  for  the  solutions  obtained  from  the  constants  of  Table  2-B  since 
both  sets  of  constants  shared  similar  solution  sensitivity  to  loading  type,  time 
increments,  and  error  tolerances.  The  default  solution  control  parameters  from 
which  these  cases  were  run  are  provided  in  Table  4-B. 


The  percent  deviation  between  the  ABAQUS  and  Mathematica  solutions 
depends  largely  on  loading  type  and  number  of  solution  increments.  Under 
isothermal  conditions,  cases  1,  2,  and  3  show  good  agreement  between  the 
ABAQUS  and  Mathematica  solutions,  as  illustrated  in  Figures  4.1,  4.2,  and  4.3, 
respectively.  The  percent  deviations  are  below  0.2%  for  all  three  cases,  which  is 
well  within  the  required  accuracy  for  most  applications. 

Table  4-A 

Verification  Test  Cases 


Test  Case 
Number 

Temp. 

(°C) 

Thermal 

Conditio 

n 

Strain 

Rate 

(lOV) 

Number 

of 

Solution 

Steps 

Solution 

File 

Appendix 

1 

25 

isotherm. 

8333 

48 

vtl.math 

3 

■  2 

650 

isotherm. 

8333 

48 

vt2.math 

4 

3 

650 

isothemx 

8333 

48 

vt3 .math 

5 

4 

25/482/25 

non- 

isotherm. 

8333 

192 

vt4.math 

6 

5 

650/760/ 

650 

non- 

isotherm 

8333 

% 

vtS.math 

7 

Table  4-B 

Default  Solution  Control  Parameters 
Parameter_ Default  Values_ Description 


INTER 

10 

number  of  initial  sub-increment 

RELAX 

0.75 

relaxation  factor 

TOLER 

10-^ 

convergence  tolerance 

XNEWDT 

1.0 

auto  time-step  size  flag 

20 


1200 


1.00 


Figure  4.1 


Figure  4.2 
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Figure  4.3  Stress-strain  response  and  solution  deviation  for  case  3. 


Two  nonisothermal  cases  are  proposed  to  evaluate  the  performance  of  the 
numerical  algorithms  tmder  combined  thermal  and  mechanical  loads.  Cases  4 
and  5  consist  of  mechanical  loads  with  two  different  temperature  ranges.  The 
temperature  range  for  case  4  is  from  25°  to  482°C,  which  covers  a  strain-rate 
insensitive  region  of  Timetal®21S.  The  temperature  profile  and  stress-strain 
response  from  case  4  are  illustrated  in  Figure  4.4,  with  the  associated  finite 
element  solution  deviation  from  Mathematica  presented  in  Figure  4.5.  Note  that 
the  error  is  minimal  during  the  isothermal  portion  of  loading,  while  the  error 
increases  slightly  with  increasing  temperature  and  then  remains  relatively 
constant. 
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The  temperature  profile  for  case  5  is  from  650°  to  760°C,  (see  Figure  4.6) 
where  the  Timetal®21S  exhibits  more  strain-rate  sensitivity  and  time  dependent 
response.  The  respective  solution  deviation,  given  in  Figure  4.7,  shows  the  finite 
element  stresses  are  well  within  0.5  percent  of  predicting  the  Mathematica  results 
for  the  isothermal  segments  of  the  loading.  The  percent  deviation  increases 
slightly  during  the  non-isothermal  periods  of  loading  and  reaches  a  single  peak 
at  2.4%.  More  solution  increments  are  required  to  improve  the  solution  accuracy, 
as  discussed  in  Section  5. 


Mechanical  Strain 


Figure  4.6  Temperature  profile  and  stress  response  for  case  5. 
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Figure  4.7  Stress-strain  response  and  solution  deviation  for  case  5. 


Solution  accuracy  and  computer  CPU  usage  for  the  test  cases  are 
siunmarized  in  Table  4-C.  For  all  the  isothermal  cases,  the  numerical  accuracy  is 
very  good  (less  than  0.5%).  The  numerical  accuracy  for  the  nonisothermal 
simulations  was  less  than  desirable.  Thus  an  investigation  was  considered  to 
optimize  numerical  accuracy  with  respect  to  computational  efficiency.  The 
findings  from  this  investigation  are  given  below. 

The  computational  speed  of  each  of  the  cases  is  quite  low  («10  CPU 
seconds).  Note  that  the  computer  usage  is  increased  when  the  simulation 
conditions  are  at  elevated  temperatures  and  slower  loading  rates.  This  increase 
in  CPU  usage  is  consistent  with  more  solution  iterations  required  when 
significant  inelastic  deformation  occurs  within  the  material. 
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I 

Table  4-C 


Summary  of  Solution  Performance 


Case 

Solution 

Increments 

Step  Size  (sec) 

CPU*  (sec) 

Maximum 
Error  (%) 

1 

48 

1.0 

3.2 

0.07 

2 

48 

1.0 

5.5 

0.14 

3 

48 

100.0 

12.6 

0.20 

4 

192 

0.5 

12.8 

1.20 

5 

96 

,  0.5 

10.0 

2.41 

^Convex,  Metaseries  2 
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SECTION  5.0 


COMPUTATIONAL  PERFORMANCE 

Case  5  from  Section  4  revealed  an  tmacceptable  level  of  solution  deviation 
between  finite  element  and  Mathematica  solutions  (2.41%).  Thus,  a  systematic 
study  was  conducted  to  determine  the  best  method  to  improve  the  accuracy  of 
the  solution.  The  performance  of  the  subroutines  was  evaluated  with  the 
tiiermomechanical  loading  of  case  5.  The  investigation  considered  three  basic 
t5^es  of  parameters  ~  number  of  solution  increments,  relaxation  factor  and 
convergence  tolerances. 

5.1  Influence  of  Number  of  Solution  Increments 


The  number  of  solution  increments  can  be  changed  at  two  different  levels 
—  major  and  minor.  The  global  convergence  of  displacements  occurs  at  the  major 
solution  increment,  NSTEP.  In  keeping  all  parameters  at  their  default  values, 
the  NSTEP  was  changed  from  96  to  192  and  384.  The  solution  accuracy  is 
improved  significantly,  from  2.41%  to  0.12  %,  as  shown  in  Figure  5.1.  For  the 
improved  accuracy,  the  CPU  usage  more  tihan  doubled,  from  10  to  23  CPU 
seconds!.  In  reference  to  Table  5-A,  cases  5, 5b  and  5c  correspond  to  NSTEP=96, 
192,  and  384,  respectively . 

At  individual  Gauss  points,  the  variable  INTER  sub-divides  the  major 
increment  into  equal  segments.  Similar  to  increasing  the  number  of  global  steps, 
increasing  the  number  of  sub-increments  from  10  to  40  improves  the  solution 
accuracy  from  2.41%  to  0.092%,  as  illustrated  in  Figure  5.2.  The  improvement  in 

!  CPU  time  on  a  Convex  Metaseries  2 
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accuracy  did  not  increase  CPU  usage  significantly.  In  reference  to  Table  5-A, 
cases  5, 5d  and  5e  correspond  to  INTER=10, 20,  and  40,  respectively. 


Mechanical  Strain 


Figure  5.1  Number  of  major  time  steps  influence  on  solution  deviation. 


Mechanical  Strain 


Figure  5.2  Number  of  sub-increments  INTER  influence  on  solution  deviation. 


28 


5.2  Influence  of  Relaxation  Factor 


In  principle,  changes  in  the  relaxation  factor  should  not  influence  the 
accuracy  of  the  solution;  however,  the  factor  should  reduce  the  number  of 
iterations  required  for  convergence.  In  the  ideal  linear  elastic  case,  a  relaxation 
factor  of  1.0  should  provide  the  quickest  convergence,  while  a  relaxation  factor  of 
0.5  is  best  suited  for  highly  non-linear  segments  of  loading.  For  case  5,  the 
relaxation  factor  does  not  significantly  improve  the  accuracy  of  the  solution,  as 
illustrated  in  Figure  5.3.  The  relaxation  factor  of  0.9  provided  a  slight 
improvement  in  CPU  usage  (9.7  seconds!  )  compared  to  that  of  RELAX=0 .75 
(10.0  seconds!  ).  Previous  experience  [20],  showed  that  a  relaxation  factor  of  0.75 
is  best  suited  for  general  usage.  In  reference  to  Table  5-A,  the  relaxation  factors 
of  0.75,  0.65  and  0.90  correspond  to  cases  5, 5f  and  5g  respectively. 


Mechanical  Strain 


Figure  5.3  Relaxation  factor  RELAX  influence  on  solution  deviation. 
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5.3  Influence  of  Solution  Tolerance 


Within  the  Bodner-Partom  iteration  algorithm,  the  tolerance  parameter 
TOLER  measures  the  change  of  stress  state  from  its  previous  guess.  A  change  of 
TOLER ,  from  10"^  to  10'^  lowers  the  level  of  error  by  a  factor  of  two,  from  2.4%  to 
1.3%,  and  with  TOLER=10'*  minimal  improvement  is  foimd,  as  illustrated  in 
Figure  5.4.  As  expected,  the  change  in  TOLER  from  lO"^  to  10'*  does  increase  the 
computer  usage  from  10.0  to  23.4  CPU  seconds^.  Values  of  TOLER  larger  than 
10'^  can  produce  non-convergent  solutions,  especially  if  the  finite  element 
configuration  contains  additional  non-linearities  (e.g.,  gap  elements).  A  smaller 
value  for  TOLER  may  be  required  to  achieve  convergence  for  these  highly  non¬ 
linear  situations.  In  reference  to  Table  5- A,  cases  5,  5h  and  5i  correspond  to 
TOLER  =  10^,  10'^  and  10“^  respectively. 


Figure  5.4  Influence  of  TOLER  on  solution  deviation. 
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Global  tolerance  GTOL  measures  the  convergence  of  displacements  and 
forces  (energy)  within  the  finite  element  mesh  for  each  major  time  step.  For  this 
investigation,  GTOL  is  set  by  the  following  input  cards: 

*CONTROLS,  PARAMETER=FIELD,  FIELD=GLOBAL 
0.0005 

Tightening  the  global  tolerance  (smaller  GTOL)  did  not  improve  the  accuracy  of 
the  solution,  in  fact  for  GTOL=5.xl0‘*,  the  solution  deviation  actually  increased 
slightly.  Under  this  tolerance,  the  CPU  usage  increased  62%,  from  10.0  to  16.2 
seconds!  when  compared  to  the  base-line  case.  In  reference  to  Table  5-A,  cases  5, 
5h  and  5i  correspond  to  GTOL  =  5x10'^  5x10^,  and  5x10’^  respectively. 
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Figure  5.5  Global  tolerance  GTOL  influence  on  solution  deviation. 
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5.4  Summary  of  Performance 


Overall,  the  results  can  be  readily  summarized  by  the  defiiution  of  the 
performance  factor: 


P  =- 


cpu 


cpu 


%D 

%D 


Eq.  10 


where  and  %D  are  the  CPU  time  and  percent  solution  deviation  for  each  test 
case.  The  values  of  and  %D’  are  the  CPU  time  and  percent  solution 

deviation  for  the  base-line  case  5.  Higher  values  of  are  indicative  of  improved 
accuracy  and/or  lower  computer  usage.  Table  5- A  presents  a  summary  of  the 
influence  of  solution  parameters  on  computer  CPU  time,  solution  deviation  and 
performance  for  the  test  cases  in  Section  5. 

Increases  in  the  number  of  solution  increments  significantly  decrease  the 
solution  deviation,  as  shown  in  Figure  5.6.  Some  improvement  is  observed  with 
decreases  in  TOLER,  with  minimal  improvement  with  tighter  global  tolerance 
GTOL.  Additional  solution  increments  and  tighter  tolerances  generally  increase 
CPU  usage,  as  illustrated  in  Figure  5.7.  The  comparison  of  overall  performance, 
P^,  as  seen  in  Figure  5.8,  revealed  that  the  best  improvement  in  accuracy  with  the 
least  debit  in  CPU  usage,  is  found  in  increases  in  INTER. 
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Table  5-A 


Summary  of  Solution  Deviation.  CPU  Usage  and  Performance  Rating 


CPU 


Case 

NSTEP 

INTER 

RELAX 

TOLER 

GTOL 

%D 

timel 

(sec) 

P. 

5 

96 

10 

0.75 

10-^ 

5x10-" 

2.41 

10.0 

1.00 

5b 

192 

10 

0.75 

10-" 

5x10'' 

1.60 

13.9 

0.98 

5c 

384 

10 

0.75 

10^ 

■saiUH 

0.077 

26.3 

0.76 

5 

96 

10 

0.75 

10^ 

5x10'" 

2.41 

5d 

96 

20 

0.75 

10-^ 

5x10-" 

1.55 

5e 

96 

40 

0.75 

lO*^ 

5x10'" 

0.092 

wsst 

Ol 

96 

10 

0.60 

10-" 

5x10'" 

10.7 

0.87 

5 

96 

10 

0.75 

10-^ 

5x10-" 

2.41 

10.0 

1.00 

% 

96 

10 

0.90 

10-" 

5x10-" 

9.70 

1.04 

5a 

96 

10 

0.75 

10-^ 

5x10'" 

2.41 

10.0 

1.00 

5h 

96 

10 

0.75 

10'" 

5x10'" 

16.2 

0.93 

5i 

96 

10 

0.75 

10-" 

5x10'" 

23.4 

0.71 

5a 

96 

10 

0.75 

10-^ 

5x10'" 

2.41 

10.0 

1.00 

5] 

96 

10 

0.75 

10^ 

5x10^ 

2.41 

11.7 

0.93 

5k 

96 

10 

0.75 

10^ 

5x10'" 

2.45 

16.2 

0.77 

■t"  CPU  time  on  a  Convex  Metaseries  2 
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CPU  Time  (sec  )  Percent  Solution  Deviation 


Case 


summary  of  solution  deviation. 


Figiure  5.7  Summary  of  CPU  usage 


SECTION  6.0 


STRUCTURE  AND  ALGORITHMS 

This  section  describes  the  structure  and  numerical  algorithms  of  the 
Bodner-Partom  subroutines  that  solve  the  constitutive  equations  of  Section  2.1 . 
Several  tables  briefly  describe  subroutine  nomenclature  and  variable  usage. 
Also,  a  description  of  the  main  subroutine's  numerical  algorithm  is  included. 
The  subroutines  described  in  this  section  are  listed  in  Appendix  A  and  are 
available  by  electromc  file  transfer. 

The  user-defined  material  subroutines  for  ABAQUS  are  written  for 
general  element  usage;  thus,  the  use  of  two  different  sets  of  subroutines  for  2- 
and  3-dimensions  is  not  required.  The  algorithms  within  ABAQUS  account  for 
the  active  number  of  normal  (NDI)  and  shear  (NSHR)  components  of  stress  and 
strains.  The  only  exception  is  for  the  plane  stress  case,  which  requires  special 
input  flag.  To  use  the  plane-stress  elements  (e.g.,  CPS8)  it  is  required  that 
IPLANE  =  1,  rather  than  I  PLANE  =  0. 
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6.1  Subroutines  and  Common  Blocks 


The  functionality  of  subroutines  found  in  usub_cab ,  f  (see  Appendix  A)  is 
summarized  in  Table  6-A.  Subroutines  UMAT  and  SDVINI  are  native  to  the 
ABAQUS  code.  Subroutines  SETPROPS ,  GETPROPS,  STGET3,  STPUT3,  and 
DBODNER  are  written  by  the  author.  The  first  four  of  these  subroutines  control 
the  values  of  material  properties  and  state  variables.  The  subroutine  DBODNER 
contains  the  numerical  algorithms  that  solve  the  Bodner-Partom  constitutive 
equations  foimd  in  Section  2.1. 


Table  6-A 

Subroutine  Names  and  Functions 


Subroutine* 

Called  From 

Function 

UMAT 

ABAQUS 

Determines  stresses  from  user-defined 
material  behavior  models 

SDVINI 

ABAQUS 

Initializes  state  variables 

SETPROPS 

UMAT 

Initializes  the  material  properties 

GETPROPS 

UMAT 

Assigns  material  properties 

STGET 

UMAT 

Retrieves  state  variables  from  working 
variable  STATEV 

STPUT 

UMAT 

Stores  state  variables  from  array 
STATEV 

DBODNER 

UMAT 

Solves  the  constitutive  equations  of  the 
Bodner-Partom  model 
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The  four  common  blocks  —  MATCONST,  BPSTATE,  MATPROP,  and 
CONTROL  —  retain  values  of  certain  variables  entering  or  leaving  the  subroutines 
in  which  they  appear.  MATCONST  transfers  material  constants  between  UMAT, 
DBODNER  and  GETPROPS.  Variables  that  make  up  common  block  MATCONST  are 
listed  in  Table  6-B.  Common  block  BPSTATE  transfers  the  values  of  the  internal 
state  variables  (e.g.,  inelastic  strain,  isotropic  hardening,  etc.)  between  all  the 
subroutines.  The  state  variables  of  common  block  BPSTATE  are  listed  in  Table  6- 
C. 


Common  block  CONTROL  transfers  the  values  of  variables  that  control 
various  aspects  of  the  solution,  e.g.,  convergence  and  time  step  size.  The  variable 
names  within  common  block  CONTROL  are  listed  Table  6-D.  Coimnon  block 
MATPROP  stores  all  the  temperature  dependent  material  properties  that  are 
initialized  in  subroutine  SETPROPS.  The  variables  contained  within  common 
block  MATPROP  are  described  in  Table  6-E. 


Table  6-B 


Internal  Variables  of  Common  Block  MATCQNST 


Variable  Name 

Symbol 

Variable  Description 

E 

E 

elastic  modulus 

AN 

n 

kinetic  (strain  rate  sensitivity)  parameter 

ZO,  Z2 

Zfl,  Zj 

fully  recovered  (minimum)  value  of  isotropic 
hardening 

Z3 

Z3 

limiting  value  of  directional  hardening 

AM2 

m. 

rate  coefficient  of  directional  hardening 

Al,  A2 

A,  A, 

thermal  recovery  coefficients  for  hardening. 

AMI 

m, 

isotropic  hardening  rate  coefficient 

Z1 

z. 

limiting  value  of  isotropic  hardening 

Rl,  R2 

thermal  recovery  exponent  for  hardening. 

DO 

Do 

limiting  inelastic  strain  rate 

ANU 

V 

Poisson's  ratio 

G 

G 

shear  modulus 

A3K 

3«K 

3  times  bulk  modulus 

DG 

dG  dT  , 
dT 

differential  shear  modulus  with  solution 
increment 

D3K 

0.  dT  , 

3  •  _  -1  dt 

dT  dt 

three  times  the  differential  bulk  modulus 
with  solution  increment 

DNU 

dv  dT 
dT  dt"^^ 

differential  Poison's  ratio  with  solution 
increment 

DM1 

5  m,  dT  , 
dT 

differential 'm,'  with  solution  increment 

Table  6-B  (continued) 

Internal  Variables  of  Common  Block  MATCONST 


Variable  Name  Symbol_ Variable  Description 


DM2 

Sm^  5  T  , 
5T 

differential 'm^'  with  solution  increment 

DRl 

dridT 

dT 

differential  rj  with  solution  increment 

DR2 

dr^dT  , 

ST  3t‘^ 

differential  r^  with  solution  increment 

DN 

d  ndT  , 
dT 

differential  'n'  with  solution  increment 

DZO 

dZodT  . 
ST  dt‘“ 

differential  with  solution  increment 

DZ2 

dZ.dT 

di 

differential  with  solution  increment 

DZ3 

dZ,dT  , 
ST  S,  * 

differential  Z3  with  solution  increment 

DAI 

d  Aid  T 

ST  St  * 

differential  Aj  with  solution  increment 

DA2 

dAidT 

ST  St  * 

differential  Aj  with  solution  increment 
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Table  6-C 


State  Variables  of  Common  Blocks  BPSTATE 


Variable  Name 

Symbol 

Variable  Description 

EIN 

e" 

ij 

components  of  inelastic  strain 

ZI 

Z. 

current  isotropic  drag  stress 

SIGEFF 

V3L 

effective  stress 

ZD 

Zo 

current  magnitude  of  directional  drag  stress 

BETA 

pi] 

components  of  directional  drag  stress 

XOUTl ,2,3 

extra  state  variables  that  can  be  used  for 
output  (especially  nice  for  debugging) 

Table  6-D 

Variables  of  Common  Block  CONTROL 

Variable  Name 

Common 

Values 

Variable  Description 

XNEWDT 

0.5 

ABAQUS  time  increment  cutting  factor 

RELAX 

0.75 

relaxation  factor  for  numerical  integration 

TOLER 

10^  to  10-" 

tolerance  for  solution  converge  for  stress  and 
state  variables 

INTER 

10 

initial  sub-time  cutting  factor 

I  PLANE 

0 

=  1  for  plane  stress  solution 

=  0  for  all  other  cases 
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Table  6-E 


Variables  of  Common  Block  MATPROP 
Variable  Name  Variable  Description 


T 

TE 

TNU 

TN 

TZ2 

TZ3 

TM2 

C5 

C6 

NTSET 


Temperature  array 

Temperature  dependent  elastic  modulus  array 

Temperature  dependent  Poisson's  ratio  array 

Temperature  dependent  material  property  n 

Temperature  dependent  material  property 

Temperature  dependent  material  property  Z3 

Temperature  dependent  material  property  m^ 

Coefficient  Cj  for  Aj  and  material  properties 
(see  Table  2) 

Exponential  Cg  for  Aj  and  A^  material  properties 
(see  Table  2-B) 

Number  of  members  in  temperature  dependent  arrays 
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6.2  Internal  Variable  Names 


This  section  summarizes  the  internal  variables  of  the  user-defined 
subroutines.  The  ABAQUS  internal  variables  that  enter  into  the  user-defined 
subroutine  UMAT  are  listed  in  Table  6-F.  The  descriptions  of  the  other  variables 
are  found  in  the  ABAQUS  user  manuals  [2].  The  primary  variables  found  in 
subroutine  DBODNER  are  summarized  in  Table  6-G. 


Table  6-F 

ABACUS  Variables  Supplied  to  UMAT  Subroutine 


Variable  Name 


STRESS 

STATEV 

DDSDDE 

STRAW 

DSTRAN 

TIME 

DTIME 

TEMP 

DTEMP 

NDI 

NSHR 

NTENS 

NSTATV 

PROPS 

COORDS 

PNEWDT 

NOEL 

NPT 

KSTEP 

KING 


Variable 

Description 

Components  of  stress 
State  variables 
Tangent  stiffness  Matrix 

Mechanical  strain  at  beginning  of  solution  increment 

Differential  mechanical  strain 

Time 

Differential  time 

Temperature 

Differential  Temperature 

Number  of  normal  components  of  stress/ strain 

Number  of  shear  components  of  stress /strain 

Total  number  of  components  of  stress/ strain 

Number  of  state  variables 

User-defined  material  properties  (see  Section  5) 

Spatial  coordinates  of  Gauss  point  location 

Flag  for  time  step  cutting 

Current  element  number 

Current  Gauss  point  niimber 

Current  ABAQUS  time  step 

Current  ABAQUS  time  increment 
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Variable  Name 

DEVEPS 
DDSDDE 
AVGEPS 
DDEPS 
DEIEFF 
DEPS 
DEPSI 
DEVSIG 
ICOUNT 
I  SUB 
RELAX 
TOLER 
SIGNEW 
SIGOLD 
SIGEST 
SIGHYD 
STRESO 
STRESS 
TIME2 
EINEST 
EINO 
BETADOT 
BETAEST 
BETAO 
U,  V 

ZIEST 


Table  6-G 

Internal  Variables  in  DBQDNER 

Variable  Description 

Deviatoric  strains 
Stiffness  matrix 

Average  (mean)  normal  strain  (bulk  strain/3) 

Sub-incremental  deviatoric  strains 

Effective  inelastic  strain  increment 

Incremental  mechanical  strain  (total  minus  thermal) 

Incremental  inelastic  deviatoric  strain 

Deviatoric  stress 

Iteration  loop  counter 

Number  of  sub-increments  within  DBODNER 
Relaxation  factor  for  new  stress  estimate 
Tolerance  for  convergence 

Effective  stress  of  previously  converged  stress  state 

Old  estimate  effective  stress 

New  effective  stress  estimate 

Hydrostatic  (mean)  stress 

Stress  of  previously  converged  stress  state 

Current  stress 

Current  time  value 

Estimated  inelastic  strains 

Inelastic  strains  of  previously  converged  stress  state 

Directional  drag  stress  rate  vector 

Estimated  directional  drag  stress  vector 

Directional  drag  stress  of  previously  converged  stress  state 

Directional  unit  vectors  for  stress  and  directional 

hardening 

Estimated  isotropic  drag  stress 
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6.3  Numerical  Algorithms 

The  algorithms  that  solve  the  Bodner-Partom  constitutive  equations 
consist  of  a  mixture  of  iteration  loops  and  sub-incremental  schemes.  This 
combination  of  iterations  and  sub-incrementation  works  well  for  the  inherently 
stiff  nature  of  the  equations  and  Timetal®21S  material  parameters.  The 
algorithms  discussed  in  this  section  are  foimd  in  subroutine  DBODNER. 

The  algorithm,  as  shown  in  Figure  6.1,  consists  of  an  iteration  loop  and  a 
sub-incremental  solution  scheme.  Prior  to  sub-incremental  integration,  state 
variable  values  assume  their  pre-incremental  values;  the  number  of  sub- 
increments  (I SUB)  initializes  to  INTER;  and  constant  rate  variables  are  scaled  by 
TFACTOR.  The  integration  of  the  sub-incremental  loop  begins  at  Step  HI.  The 
primary  iteration  loop  starts  at  step  IH.D  and  converges  on  stress,  inelastic  strain, 
isotropic  hardening  and  the  directional  hardening  parameters.  Non-convergent 
solutions,  as  defined  by  a  maximum  limit  set  on  iteration  steps,  return  to  the 
beginning  of  the  sub-incremental  integration  (Step  I)  with  an  increase  in  the 
number  of  sub-increments,  I  SUB.  When  the  maximiim  munber  of  sub- 
incremental  step  equals  128,  the  solution  solver  ceases  operation  and  prints  a 
diagnostic  debugging  output. 


I.  Store  variable  at  the  beginning  of  solution  increment 

A.  state  variables 

B.  temperature-dependent  material  parameters 

n.  Initialize  variables  for  sub-incremental  solution  step  cutting  of  I  SUB 

A.  initialize  estimated  new  values  for  state  variables 

B.  restore  material  parameters  from  Step  I 

C.  determine  new  material  parameters  rates  by  factors  of  I  SUB 

D.  restore  state  variables  to  values  from  Step  I 

m  Begin  sub-incremental  solution  integration 

A.  update  aU  temperature-dependent  parameters  to  end  of  sub- 
incremental  step 

B.  update  deviatoric  and  mean  mechanical  strains 

C.  step  iteration  coimter  ICOUNT 

D.  begin  primary  iteration  loop 

1.  estimate  inelastic  strain  increment  (engineering 
shear  strains  are  computed) 

2.  compute  new  stress  state  and  inelastic  work  rate 

3.  compute  new  estimate  for  beta 

4.  compute  new  estimate  for  isotropic  hardening 

5.  check  for  convergence  of  effective  stress,  incremental 
inelastic  strains,  isotropic  and  directional  hardening 

a.  if  converged,  then  continue  witin  conclude  current  sub- 
increment  via.  STEP  ni.E 

b.  if  not,  make  new  estimate  for  effective  stress  and 
increase  iteration  count  ICOUNT  by  1. 

c.  check  for  excessive  iteration  coimt 

(1)  if  excessive,  then  increase  cutting  factor  I  SUB 
by  factor  of  2.0  and  precede  to  Step  H.  A 

(2)  if  not,  to  Step  IH.D 

E.  Update  converged  solution  with  estimates,  then  continue  at 
STEP  in.A 

IV.  Complete  aU  sub-incremental  cycles  and  then  return  to  UMAT 

Figure  6.1  Numerical  algorithm  for  solving  of  Bodner-Partom  equations. 
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APPENDIX  A 


Listing  Of  User-Defined  Subroutine 

*USER  SUBROUTINE 

SUBROUTINE  UMAT ( STRESS , STATEV, DDSDDE , SSE , SPD , SCO , 

1  RPL , DDSDDT , DRPLDE , DRPLDT , STRAN , DSTRAN , 

2  TIME, DTIME, TEMP, DTEMP,PREDEP,DPRED,MATERL,NDI,NSHR,NTENS, 

3  NSTATV , PROPS , NPROPS , COORDS , DROT , PNEWDT , CELENT , 

4  DFGRDO , DFGRDl , NOEL , NPT , KSLAY , KSPT , KSTEP , KING ) 

C 

INCLUDE  ' ABA_PARAM . INC ' 

CHARACTER* 8  CMNAME 
C 

DIMENSION  STRESS (NTENS) ,STATEV{NSTATV) , 

1  DDSDDE ( NTENS , NTENS ) , DDSDDT { NTENS ) , DRPLDE ( NTENS ) , 

2  STRAN(NTENS) , DSTRAN (NTENS) ,TIME{2) ,PREDEF{1) ,DPRED{1) , 

3  PROPS (NPROPS) ,COORDS(3) ,DROT(3,3)  , 

4  DFGRDO (3, 3), DFGRDl (3, 3) 

C 

PARAMETER  (ONE=l . ODO , TWO=2 . ODO , THREE=3 . ODO , SIX=6 . ODO ) 

C 

C  THIS  SOFTWARE  AND  ANY  ACCOMPANYING  DOCUMENTATION  IS  RELEASED  AS 

C  IS.  THE  U.S.  GOVERNMENT,  ITS  CONTRACTORS  AND  THEIR  SUBCONTRACTORS 

C  MAKE  NO  WARRANTY  OF  ANY  KIND,  EXPRESS  OR  IMPLIED,  CONCERNING  THIS 
C  SOFTWARE  AND  ANY  ACCOMPANYING  DOCUMENTATION,  INCLUDING,  WITHOUT 
C  LIMITATION,  ANY  WARRANTIES  OF  MERCHANTABILITY  AND  FITNESS  FOR  A 
C  PARTICULAR  PURPOSE.  IN  NO  EVENT  WILL  THE  U.S.  GOVERNMENT,  ITS 
C  CONTRACTORS  AND  THEIR  SUBCONTRACTORS  BE  LIABLE  FOR  ANY  DAMAGES, 

C  INCLUDING  LOST  PROFITS,  LOST  SAVINGS  OR  OTHER  INCIDENTAL  OR 
C  CONSEQUENTIAL  DAMAGES  ARISING  OUT  OF  THE  USE,  OR  INABILITY  TO  USE 

C  THIS  SOFTWARE  OR  ANY  ACCOMPANYING  DOCUMENTATION,  EVEN  IF  INFORMED 

C  IN  ADVANCE  OF  THE  POSSIBILITY  OF  SUCH  DAMAGES. 

C 

C  "This  software  is  being  used  at  the  user's  own  risk:  Neither  the 
C  Government  Agency  nor  its  contractors  assure  software ' s  accuracy 
C  or  its  appropriate  use . " 

C 

c  Directional  Bodner-Partom  Material  Constants: 
c 

c  Temperature  Independent  Constants 

c 

c  do,  Z1 

c 

c  Temperature  Dependent  Constants 

c 

c  an,  ami,  am2,  zO,  z2,  z3,  al,  a2 ,  rl,  r2 

c 

c  Temperature  Differentials 

c 

c  dn,  dml,  dm2,  dzO,  dz2,  dz3,  dal,  da2,  drl,  dr2 

c 

COMMON  /matconst/  E,  G,  a3k,  anu,  dg,  d3k,  dnu,  dO, 

&  an,  ami,  am2,  zO,  zl,  z2,  z3,  al,  a2,  rl,  r2. 

Sc  dn,  dml,  dm2,dz0,dz2,dz3,dal,da2,drl,dr2 

C 

c  common  block  for  state  varibles 
c 

common  /bpstat/  ein( 6) , zi, sigef f ,betaef f , epef f , 

&  zd,  beta (6),  xoutl,  xout2,  xout3 
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c 

common  /control/  xnewdt.  relax,  toier,  inter,  Iplane 
c 

c 


c 

Set 

Iplane  = 

0  for 

Three- D 

c 

for 

Plane  Strain 

c 

for 

Ax i symmetric 

c 

Set 

Iplane  - 

1  for 

Plane  Stress 

^  ****★★****★********•*************♦-********•*************** 

c 

c 

if(inter  .eq.  0)  then 
call  setprops ( PROPS , NPROPS ) 
endi  f 


c 

c  recall  state  varibles 
c 

call  stget ( STATEV, NSTATV; 
c 

call  getprops ( TEMP , DTEMP ) 
c 

time2  =  time ( 1 ) 
c 

call  dbodner ( stran, ds tran , stress , ddsdde, dt ime , time2 , 
$  pnewdt ,  ndi , nshr , ntens , npt ,  kstep) 

c 


c 

c  If  Bodner  Partom  solution  did  not  converge 
c  then  send  a  cut  the  major  t ime  step  defined  by 

c  abaqus 

c 

call  stput (STATEV, NSTATV) 


Elastic- Plastic  flow  correction  for  DDSDDE 


c 

(Currently  Di 

sabied  -  see  Joe,  tor  Details) 

c 

c 

sigave  =  0 

.0 

c 

do  77  i  ^ 

1,  NDI 

c 

sigave  =  sigave  +  stress(i) 

c  77 

continue 

c 

sigave  =  s 

igave  /  THREE 

cc 

c 

do  78  i  = 

1 ,  NTENS 

c 

i f  ( sigef f 

•  It.  1 . Oe  - 10 )  then 

c 

f 1 ow ( I ) 

=  0.0 

c 

else 

c 

if(i  .gt.  ndi)  then 

c 

f 1 ow ( i ) 

=  s tress ( i )/ s igef f 

c 

else 

c 

flow ( i ) 

=  ( St ress ( 1 )  sigave) /sigef f 

c 

endi  f 

c 

endi  f 

c  78 

continue 

cc 

c 

YIELD 

SIGEFF 

cc 

c 

EFFG3 

=  THREE*G 
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EFFHRD 


EFFG3 


c 

cc 

c  DO  250  K1=1,NTENS 

c  DO  240  K2=1,NTENS 

C  DDSDDE (K2, K1 ) =DDSDDE( K2 , K1 ) +FLOW (K2 ) *FLOW{Kl ) 

c  1  * (EFFHRD-EFFG3 ) 

c  240  CONTINUE 

c  250  CONTINUE 

cC 

RETURN 

END 

C 

SUBROUTINE  SDVINI ( STATEV , COORDS , NSTATV , NCRDS , NOEL , NRT , LAYER, KSPT) 
C 

INCLUDE  ' ABA_PARAM . INC 
C 

DIMENSION  STATEV ( NSTATV ) ,  COORDS ( NCRDS ) 
c 

COMMON  /matconst/  E,  G.  a.iV,  anu,  dg,  d3k,  dnu ,  dO , 

&  an,  ami,  am2 ,  zO,  zl,  22,  zJ,  al.  a2 ,  rl ,  r2, 

&  dn,  dml ,  dm2 , dzO . dz2 , dz3 , dal , da2 , dr 1 , dr2 

c 

common  /control/  etthrd,  xnewdt ,  relax,  Coier,  inter,  Iplane 
c 

inter  =  0 
c 

statev{l)  -  0.0 
statev{2)  =  0.0 
statev(3)  0.0 
statev(4)  =  0.0 
statev(5)  =  0.0 
statev(6)  =  0.0 
statev(7)  -  0.0 
statev ( 8 )  =  0.0 

statevO)  =  0.0 
statev ( 10 )  -  0.0 

statev(ll)  =  0.0 
statev(12)  =  0.0 
c 

c  Initializing  ZI  to  ZO 
c 

c  tor  Neu's  Constants,  Table  2-B  [4] 
c 

c25 


c25 

statev (13) 

=  1550 

c 

c650 

c6b0 

statev (13) 

=  1000 

c 

c900 

statev (13 ) 

=  300 

c 

c  for  Sander's  Constants,  Table  2-C  [18] 


c 

c25 


c25 

c650 

statev (13) 

-  2950 

c650 

c815 

statev ( 13 ) 

=  200 

c815 

c 

statev ( 13 ) 

-  150 

,52 


o  n  o 


statev(14)  =  0.0 
statev(15)  =  0.0 
statev(16)  =  0.0 
statev(17)  =  0.0 
c 

return 
end 
C 

subroutine  dbodner ( eps , deps , stress , ddsdde , dtime , tinie2 , 

$  pnewdt,  ndi,nshr,ntens,npt,  kstep) 

c 

Q  *************************************************** 

c 
c 
c 
c 
c 
c 
c 
c 
c 
c 
c 
c 
c 
c 

Q  *************************************************** 

c 

C  THIS  SOFTWARE  AND  ANY  ACCOMPANYING  DOCUMENTATION  IS  RELEASED  AS 

C  IS.  THE  U.S.  GOVERNMENT,  ITS  CONTRACTORS  AND  THEIR  SUBCONTRACTORS 

C  MAKE  NO  WARRANTY  OF  ANY  KIND,  EXPRESS  OR  IMPLIED,  CONCERNING  THIS 
C  SOFTWARE  AND  ANY  ACCOMPANYING  DOCUMENTATION,  INCLUDING,  WITHOUT 
C  LIMITATION,  ANY  WARRANTIES  OF  MERCHANTABILITY  AND  FITNESS  FOR  A 
C  PARTICULAR  PURPOSE.  IN  NO  EVENT  WILL  THE  U.S.  GOVERNMENT,  ITS 
C  CONTRACTORS  AND  THEIR  SUBCONTRACTORS  BE  LIABLE  FOR  ANY  DAMAGES, 

C  INCLUDING  LOST  PROFITS,  LOST  SAVINGS  OR  OTHER  INCIDENTAL  OR 
C  CONSEQUENTIAL  DAMAGES  ARISING  OUT  OF  THE  USE,  OR  INABILITY  TO  USE 

C  THIS  SOFTWARE  OR  ANY  ACCOMPANYING  DOCUMENTATION,  EVEN  IF  INFORMED 

C  IN  ADVANCE  OF  THE  POSSIBILITY  OF  SUCH  DAMAGES. 

C 

C  "This  software  is  being  used  at  the  user's  own  risk:  Neither  the 
C  Government  Agency  nor  its  contractors  assure  software's  accuracy 
C  or  its  appropriate  use." 

C 

INCLUDE  'ABA_PARAM.INC' 

C - 

c 

real*8  Deps{ntens),  STRESS { ntens ) ,  DEPSI(ntens) 

real*8  eps (ntens),  ddeps (ntens) ,  DEVSIG (ntens) 

real*8  STRESO (ntens) ,  einest (ntens) ,  einO (ntens) 

real*8  betaest (ntens) ,  betaO(ntens) 

real*8  epsO (ntens),  ddsdde (ntens, ntens) 


COMMON  /matconst/  E,  G,  a3k,  anu,  dg,  d3k,  dnu,  dO, 
&  an,  ami,  am2,  zO,  zl,  z2,  z3,  al,  a2,  rl,  r2, 

&  dn,  dml,  dm2,dz0,dz2,dz3,dal,da2,drl,dr2 

c 

common  /bpstat/  ein(6),  zi,  sigeff,  betaeff,  epeff, 

&  zd,  beta (6),  xoutl,  xout2,  xout3 


Copyright  1995  The  University  of  Dayton 
All  rights  reserved. 

This  material  may  be  reproduced  by  or  for  the 
U.S.  Govenment  pursuant  to  the  copyright  license 
under  the  clause  at  DEARS  252.227-7013  (Oct  1988) 

Programmed  by: 

Joseph  L.  Kroupa  (513)  255-9871 
University  of  Dayton  Research  Institute 
Dayton,  Ohio  45469-0128 
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coitmvon  /control/  effhar,  xnewdt,  relax,  toler,  inter,  Iplane 


xsqrt3  =  sqrt(3.) 
store  old  state  varibles 


DO  10  I=l,ntens 
stresO(I)  =  stress (I) 
einO(I)  =  ein(I) 
beta0{i)  =  beta(i) 
epsO(i)  =  eps(i) 
CONTINUE 


do  11  kl=l,NTENS 
do  11  k2=l,NTENS 
DDSDDE(k2,kl)  =  0.0 
continue 

ziO  =  zi 

zdO  =  zd 

store  old  material  constants 

gO  =  g 

a3k0  =  a3k 
anuO  =  anu 
zOO  =  zO 
z20  =  z2 

z30  =  z3 

alO  =  al 
a20  =  a2 

ami  0  =  ami 
am20  =  am2 
rlO  =  rl 
r20  *  r2 

initialize  other  varibles 

isub  =  inter 
i submax  =  inter  *  64 
idbug  =  0 

isecond  =  0 


initialize  variables  for  sub-time  cutting 


continue 

sigest  =  sigeff 

zi  =  ziO 

ziest  =  ziO 

zold  =  ziO  +  zdO 

zd  =  zdO 

restore  old  material  constants 
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z3 

= 

z30 

al 

s: 

alO 

a2 

= 

a20 

ami 

amlO 

am2 

= 

am20 

rl 

rlO 

r2 

= 

r20 

time 

= 

time  2 

c 

c  update  new  rate  material  constants  with  subincrement  isub 
c 

tfactor  =  1./  dfloat(isub) 
c 


dtsub 

dtime 

★ 

tfactor 

dgsub 

= 

dg 

it 

tfactor 

d3ksub 

= 

d3k 

it 

tfactor 

dnusub 

= 

dnu 

it 

tfactor 

dnsub 

= 

dn 

it 

tfactor 

dzOsub 

= 

dzO 

it 

tfactor 

dz2sub 

= 

dz2 

* 

tfactor 

dz3sub 

= 

dz3 

★ 

tfactor 

dal sub 

= 

dal 

★ 

tfactor 

da2sub 

= 

da2 

it 

tfactor 

dmlsub 

s= 

dml 

it 

tfactor 

dm2 sub 

dm2 

it 

tfactor 

dr Isub 

= 

drl 

it 

tfactor 

dr  2 sub 

dr2 

it 

tfactor 

c 

c  compute  devitoric  stress 
c 

sighyd  =0. 

do  19  i  =  1,  ndi 

sighyd  =  sighyd  +  stresO(i) 

19  continue 

sighyd  =  sighyd  /  3.0 
c 

do  2  0 ,  i  =  1 , ntens 
c 

ein(I)  =  einO(i) 

beta{i)  =  betaO(i) 

betaest(i)  =  betaO(i) 

eps (i)  =  epsO (i) 

ddeps{i)  =  deps(i)  *  tfactor 
c 

IF  (I.gt.ndi)  THEN 
f actl  =  0 . 

ELSE 

factl  =  1. 

ENDIF 

c 

devsig(I)  =  stresO(I)  -  factl*sighyd 
c 

20  CONTINUE 
c 

c 

DO  200  Jsub=l,isub 
c 

c  update  all  temperature  dependent  material  constants  to 
c  end  of  subtime  increment  step 
c 

g  -  g  +  dgsub 
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c 


c 


a3k 

= 

a3k 

anu 

= 

anu 

an 

= 

an 

zO 

= 

zO 

z2 

= 

z2 

z3 

= 

z3 

al 

= 

al 

a2 

= 

a2 

ami 

= 

ami 

am2 

= 

am2 

rl 

= 

rl 

r2 

= 

r2 

time 

= 

time 

+  dSksub 
+  dnusub 
+  dnsub 
+  dzOsub 
+  dz2sub 
+  dzSsub 
+  dal sub 
+  da2sub 
+  dmlsub 
+  din2  sub 
+  dr 1 sub 
+  dr2sub 
+  dtsub 


EFFG2  =  2.0*G 

EFFG3  =  3.0*G 


if(Iplane  .eq.  1)  then 

EFFLAM  =  EFFG2  *  (anu/ (1 . -anu) ) 

else 

EFFLAM  =  (A3K-EFFG2) /3.0 
endif 


DO  220  K1=1,NDI 
DO  210  K2=1,NDI 

DDSDDE (K2 , K1 ) =EFFLAM 
210  CONTINUE 

DDSDDE (K1 , K1 ) =EFFG2+EFFLAM 
220  CONTINUE 

DO  230  K1=NDI+1,NTENS 
DDSDDE (K1,K1)=G 
2,30  CONTINUE 

eff_el  =  0 
do  30  i  =  l,ntens 
if(i  .gt.  ndi)  then 

eff_el  =  eff_el  +  ddeps (i) *ddeps (i) *2 
else 

eff_el  =  eff_el  +  ddeps (i) *ddeps (i) 
endi  f 

eps(i)  =  eps(i)  +  ddeps (i) 
einest(I)  =  ein(I) 

30  continue 

eff_el  =  sqrt{0.5*eff_el) 

icount  =  0 
deiold  =  0 . 


c 

ztot  =  ziest  +  zd 

c 

300 

continue 

c 

c  For 
c  Set 

Elastic  Test  case 
plasticity  off 

c 

cof  f 

dO  =  0.0 

c 

IF(sigest  .It.  l.OE-30)  THEN 
sigest  =  l.Oe-30 
deieff  =  0.0 
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xlam  =  0.0 
ELSE 

xtmpl  =  (ztot/sigest) **2 
xtmp2  =  -0 . 5*xtmpl**an 
deieff  =  d0*exp (xtmp2 ) 
xlam  =  xsqrt3  *  deieff /sigest 
ENDIF 
c 

deieff  =  deieff  *  dtsub 
c 

effhrd  =  0 

if  {eff_el  .gt.  0.0)  effhrd  =  deief f /ef f_el 

xoutl  =  effhrd 
c 

DO  40  I=l,ntens 

IF  (I.gt.ndi)  THEN 
factl  =  2 
fact2  =  0 
ELSE 

factl  =  1 
fact2  =  1 
ENDIF 
C 

c  estimate  plastic  strains  and  stresses 
c  (Engineering  Plastic  Shear  Strains  are  Computed) 
c 

DEPSI(I)  =  xlam  *  devsig(I)  *  factl 
c 

c  Defined  thermal  differentail  terms 
c 

theta3  =  0.0 
c 

einest(i)  =  ein(i)  +  (depsi (i) *dtsub)  +theta3 
c 

40  CONTINUE 
c 
c 

ssum  =  0. 
pwdot  =  0 . 
sigavg  =  0. 
c 

DO  41  I=l,ntens 

IF  (I.gt.ndi)  THEN 
factl  =  2 
fact2  =  0 
ELSE 

factl  =  1 
fact2  =  1 
ENDIF 
C 

stress  (I)  =  0.0 

do  39  j  =  1,  NTENS 

stress(I)  =  stress(I)  +  ddsdde (i, j ) * (eps ( j ) -einest ( j ) ) 
39  continue 
c 

pwdot  =  pwdot  +  (stress (I) *depsi (I) ) 
ssum  =  ssum  +  factl  *  stress (I)**2 

sigavg  =  sigavg  +  fact2  *  stress (I) 

c 

c  Extra  output,  nice  for  Debugging 
c 
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if(i  .eq.  2)  then 
xoutl  =  depsi(I) 
xout2  =  stress (I) 
endif 

c 

41 

continue 

c 

ssum  =  sqrt (ssum) 

sigavg  =  sigavg/3.0 

c 

bnew  =  0.0 
signew  =  0.0 

c 

DO  42  I=l,ntens 

IF  (I.gt.ndi)  THEN 
factl  =  2 
fact2  =  0 
ELSE 

factl  =  1 
fact2  =  1 
ENDIF 
c 

bnew  =  bnew  +  factl  *  beta(i)**2 

devsig(i)  =  stress (i)  -  fact2  *  sigavg 
signew  =  signew  +  factl  *  devsig(I)**2 
c 

42  continue 
c 

bnew  =  sqrt(bnew) 
c 

c  For  plane  stress  solution 

c  Add  extra  deviatoric  stress  for  third  plane 
c 

ifdplane  .eq.  1)  signew  =  signew  +  sigavg**2 
c 

signew  =  sqrt (1 . 5*signew) 
c 

zd  =  0 . 

do  44  i  =  l,ntens 
IF  (I.gt.ndi)  THEN 
fact=2 
ELSE 
fact=l 
ENDIF 
c 

c  compute  drag  stress  vectors 
c 

V  =  beta(i) 

if (bnew  .ge.  l.E-30)  v  =  beta(i)/  bnew 
u  =  stress (i) 

if(ss\im  .ge.  l.E-30)  u  =  stress(i)  /  ssum 
c 

bterml  =  ((z3  *  u)  -  betaest(i))  *  pwdot 
bterm2  =  zl  *  v  * ( (bnew/zl) **r2) 
theta2  =  betaest(i)  *  dz3sub  /  z3 

betadot  =  am2*bterml  -  a2*bterm2 
betaest(i)  =  beta(i)  +  betadot*dtsub  +  theta2 
c  write (6,*)  'flag  6' 

c 

zd  =  zd  +  betaest(i)  *  u  *  fact 
c  write (6,*)  'flag  7' 
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o  o  n  non 


44 

continue 

c 

46 

zterml  =  zl* (abs (ziest-z2 ] 

1  /zl) **rl 

c 

thetal  =  ( (zi-z2) / (zl-z2) ) 

*  dz2sub 

c 

ZiDOT  =  ami* (zl-ziest) *pwdot 
&  -  al*zterml 


COMPUTE  UPDATED  Ztot  VALUE 

ziest  =  zi  +  zidot*dtsub  +  thetal 
ztot  =  zd  +  ziest 

if (ztot  .It.  toler)  ztot=  toler 

NOW  INVESTIGATE  CONVERGENCE 
c 

serror  =  dABS (sigest-signew) 
sigmax  =  dmaxl ( signew, sigest) 

IF(sigmax  .gt.  toler)  serror  =  serror  /  sigmax 
c 

deperr  =  dabs (deief f-deiold) 
depmax  =  dmaxl (deief f, deiold) 

IF (  depmax  .gt.  toler)  deperr  =  deperr  /  depmax 
c 

zerror  =  dABS (ztot-zold) 
zmax  =  dmaxl (ztot, zold) 

IF(  zmax  .gt.  toler)  zerror  =  zerror  /  zmax 
c 

error  =  serror  +  deperr  +  zerror 
c 

c  - 

c 

if(idbug  .eq.  1)  then 

write(6,*)  ' icount  isub  errors  s  ep  ztot  ',icount, 
write (6,*)  sigest,  zd,  serror , zerror , deperr 
endif 
c 

c - 

c 

c  check  for  convergence 
c 

IF(  error  .gt.  toler)  THEN 
c 

c  nonconvergent  solution  arrives  here 
c 

sigest  =  sigest  *  (1. -relax)  +  signew  *  relax 
c  sigest  =  signew 

deiold  =  deieff 
zold  =  ztot 
c 

c  update  and  check  convergence  count 
c 

icount  =  icount+1 
c 

IF  (icount  .ge.  10  .or.  error  .ge.  2.9)  THEN 
c 

c  restart  iteration  scheme  with  new  time-step  cut 
c  make  additional  time-step  cuts 
c 


isub 
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isub  =  isub  *  2 


terminate  nonconverge  solutions 
if  (isub  .gt.  i submax)  then 


write (6, *) 

write ( 6 , *) 

'  Bodner 

Solution 

write ( 6 , * ) 

write { 6 , * ) 

'  icount 

jsub  isub 

write ( 6 , * ) 

'  error 

,  error 

write ( 6 , * ) 

'  depsi 

.  depsi 

write ( 6 , *) 

'  ein  ' , 

ein 

write ( 6 , *) 

'  einest 

■ ,  einest 

write (6, *) 

'  dtsub 

,  dtsub 

write ( 6 , * ) 

' devsig 

,  devsig 

write (6 , * ) 

' stress 

,  stress 

write (6, *) 

'  zd  ■  , 

zd 

write (6, *) 

'  ztot  ' 

,  ztot 

stop 

icount,  jsub,  isub 


endif 

restart  iteration  scheme  with  new  time-step  cut 
goto  100  . 
endif 
goto  300 
■else 

convergence  solution  arrives  here 

zi  =  ziest 

do  60  i  =  1,  ntens 
ein(I)  =  einestd) 


beta (I) 
continue 


=  betaest(I) 


endif 

end  of  subtime  loops 

)  continue 

epeff  =  0.0 

do  70  i  =  1,  ntens 

IF  (I.gt.ndi)  THEN 
f act2  =  2 . 

ELSE 

fact2  =  1'. 

ENDIF 

epeff  =  epeff  +  ein{i)  *  ein(i)  *  fact2 
continue 

epeff  =  sqrt (2 . 0*epef f /3 . 0) 

SIGEFF  =  sigest 
betaeff  =  bnew 


c  adjust  increment  size  parameter 
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c 

pnewdt  =1.0 
c 

if  (isub  .eq.  inter  .or.  effhrd  .It.  0.80)  then 
pnewdt  =  xnewdt 
endif 
c 

RETURN 

END 

c 

subroutine  stget (Array, Nstate) 
c 

Q  *****************  ***********************  *********** 

C 

c  Copyright  1995  The  University  of  Dayton 

c  All  rights  reserved, 

c 

c  This  material  may  be  reproduced  by  or  for  the 
c  U.S.  Govenment  pursuant  to  the  copyright  license 
c  under  the  clause  at  DEARS  252.227-7013  (Oct  1988) 
c 

c  Programmed  by: 
c 

c  Joseph  L.  Kroupa  (513)  255-9871 
c  University  of  Dayton  Research  Institute 
c  Dayton,  Ohio  45469-0128 
c 

Q  *★■*★*'***★★*★★★**★**★★*******★******★**★***★★****★** 

c 

C  THIS  SOFTWARE  AND  ANY  ACCOMPANYING  DOCUMENTATION  IS  RELEASED  AS 

C  IS.  THE  U.S.  GOVERNMENT,  ITS  CONTRACTORS  AND  THEIR  SUBCONTRACTORS 

C  MAKE  NO  WARRANTY  OF  ANY  KIND,  EXPRESS  OR  IMPLIED,  CONCERNING  THIS 
C  SOFTWARE  AND  ANY  ACCOMPANYING  DOCUMENTATION,  INCLUDING,  WITHOUT 
C  LIMITATION,  ANY  WARRANTIES  OF  MERCHANTABILITY  AND  FITNESS  FOR  A 
C  PARTICULAR  PURPOSE.  IN  NO  EVENT  WILL  THE  U.S.  GOVERNMENT,  ITS 
C  CONTRACTORS  AND  THEIR  SUBCONTRACTORS  BE  LIABLE  FOR  ANY  DAMAGES, 

C  INCLUDING  LOST  PROFITS,  LOST  SAVINGS  OR  OTHER  INCIDENTAL  OR 
C  CONSEQUENTIAL  DAMAGES  ARISING  OUT  OF  THE  USE,  OR  INABILITY  TO  USE 

C  THIS  SOFTWARE  OR  ANY  ACCOMPANYING  DOCUMENTATION,  EVEN  IF  INFORMED 

C  IN  ADVANCE  OF  THE  POSSIBILITY  OF  SUCH  DAMAGES. 

C 

C  "This  software  is  being  used  at  the  user’s  own  risk:  Neither  the 
C  Government  Agency  nor  its  contractors  assure  software ' s  accuracy 
C  or  its  appropriate  use." 

C 


c 

c 


c 


INCLUDE  ■ABA_PARAM.INC’ 
real*8  array (Nstate) 


common  /bpstat/  ein(6) ,  zi,  sigeff,  betaeff ,  epeff , 
&  zd,  beta (6),  xoutl,  xout2,  xout3 


ein (1) 
ein(2) 
ein (3) 
ein (4) 
ein(5) 
ein( 6) 
beta ( 1) 
beta (2) 
beta { 3 ) 


array ( 1 ) 
array ( 2 ) 
array ( 3 ) 
array ( 4 ) 
array ( 5 ) 
array ( 6 ) 
array ( 7 ) 
array ( 8 ) 
array ( 9 ) 
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beta (4) 
beta (5) 
beta (6) 
zi 

sigef f 
epef  f 
betaef f 
zd 
c 

return 

end 

c 

subroutine  stput (Array, Ns tate) 
c 
c 

Q  *  **  *  *********************************************** 

c 

c  Copyright  1995  The  University  of  Dayton 

c  All  rights  reserved, 

c 

c  This  material  may  be  reproduced  by  or  for  the 
c  U.S.  Govenment  pursuant  to  the  copyright  license 
c  under  the  clause  at  DFARS  252.227-7013  (Oct  1988) 
c 

c  Programmed  by : 
c 

c  Joseph  L.  Kroupa  (513)  255-9871 
c  University  of  Dayton  Research  Institute 
c  Dayton,  Ohio  45469-0128 
c 

c 

C  THIS  SOFTWARE  AND  ANY  ACCOMPANYING  DOCUMENTATION  IS  RELEASED  AS 

C  IS.  THE  U.S.  GOVERNMENT,  ITS  CONTRACTORS  AND  THEIR  SUBCONTRACTORS 

C  MAKE  NO  WARRANTY  OF  ANY  KIND,  EXPRESS  OR  IMPLIED,  CONCERNING  THIS 
C  SOFTWARE  AND  ANY  ACCOMPANYING  DOCUMENTATION,  INCLUDING,  WITHOUT 
C  LIMITATION,  ANY  WARRANTIES  OF  MERCHANTABILITY  AND  FITNESS  FOR  A 
C  PARTICULAR  PURPOSE.  IN  NO  EVENT  WILL  THE  U.S.  GOVERNMENT,  ITS 
C  CONTRACTORS  AND  THEIR  SUBCONTRACTORS  BE  LIABLE  FOR  ANY  DAMAGES, 

C  INCLUDING  LOST  PROFITS,  LOST  SAVINGS  OR  OTHER  INCIDENTAL  OR 
C  CONSEQUENTIAL  DAMAGES  ARISING  OUT  OF  THE  USE,  OR  INABILITY  TO  USE 

C  THIS  SOFTWARE  OR  ANY  ACCOMPANYING  DOCUMENTATION,  EVEN  IF  INFORMED 

C  IN  ADVANCE  OF  THE  POSSIBILITY  OF  SUCH  DAMAGES. 

C 

C  "This  software  is  being  used  at  the  user's  own  risk:  Neither  the 
C  Government  Agency  nor  its  contractors  assure  software's  accuracy 
C  Or  its  appropriate  use." 

C 


=  array (10) 
=  array (11) 
=  array (12) 
=  array (13) 
=  array (14) 
=  array (15) 
=  array (16) 
=  array (17) 


c 

c 


INCLUDE  'ABA_P ARAM. INC 


real*8  array (Nstate) 


common  /bpstat/ 


zd. 

array ( 1 ) 

=  ein(l) 

array ( 2 ) 

=  ein(2) 

array ( 3 ) 

=  ein(3) 

array ( 4 ) 

=  ein(4) 

array ( 5 ) 

=  ein(5) 

array ( 6 ) 

=  ein(6) 

ein( 6) ,  zi, 
beta (6)  , 


sigeff,  betaeff,  epeff, 
xoutl ,  xout2 ,  xout3 
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array ( 7 ) 

= 

beta ( 1) 

array ( 8 ) 

= 

beta (2 ) 

array ( 9 ) 

= 

beta (3 ) 

array (10) 

= 

beta (4) 

array (11) 

= 

beta ( 5 ) 

array ( 12 ) 

= 

beta (6) 

array (13 ) 

= 

zi 

array  (14)' 

= 

sigef f 

array (15) 

=: 

epef  f 

array (16) 

= 

betaef f 

array (17) 

= 

zd 

array (18) 

= 

xoutl 

array (19) 

= 

xout2 

array (20) 

xout3 

return 

end 

c 

subroutine  setprops ( PROPS , NPROPS ) 
c 
c 

c 

c  Copyright  1995  The  University  of  Dayton 

c  All  rights  reserved, 

c 

c  This  material  may  be  reproduced  by  or  for  the 
c  U.S.  Govenment  pursuant  to  the  copyright  license 
c  under  the  clause  at  DFARS  252.227-7013  (Oct  1988) 
c 

c  Programmed  by : 
c 

c  Joseph  L.  Kroupa  (513)  255-9871 
c  University  of  Dayton  Research  Institute 
c  Dayton,  Ohio  45469-0128 
c 

c 

INCLUDE  ' ABA_PARAM.INC' 
c 

PARAMETER  (MAXP=20) 

common  /matprop/  T(MAXP) ,  TE(MAXP) ,  TNU(MAXP) ,  TN(MAXP) , 

&  TZ2 (MAXP) ,  TZ3 (MAXP) ,  TMl (MAXP) ,  TM2 (MAXP) , 

&  TRl(MAXP),  TR2{MAXP),  TINTER  (MAXP), 

&  C51,  C61,  C52,  C62,  NTSET 

c 

COMMON  /matconst/  E,  G,  a3k,  anu,  dg,  d3k,  dnu,  dO , 

&.  an,  ami,  am2,  zO,  zl,  z2,  z3,  al,  a2,  rl,  r2 , 

&  dn,  dml,  dm2 ,  dzO ,  dz2 ,  dz3  ,  dal ,  da2 ,  drl ,  dr2 

c 

common  /control/  effhrd,  xnewdt,  relax,  toler,  inter,  Iplane 
c 

real*8  Props (NPROPS) 
c 
c 

NOUT=6 

c 

c  Extract  Control  Variables 

c 

Iplane  =  IFIX ( PROPS ( 1 ) ) 

NTIP  =  IFIX(PROPS(2) ) 
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NTSET  =  IFIX(PR0PS(3) ) 
c 

write(NOUT,  '{//.34h  Input  Control  Variables  ,11)') 

write (NOUT, ' (15h  IPLANE  ,  2x,  15  )')  iplane 

write (NOUT, ' (15h  NTIP  ,  2x,  15  ) ' )  NTIP 

write (NOUT, ' (15h  NTSET  ,  2x,  15  )')  NTSET 

c 

do  =  Props (4) 
zl  =  Props (5) 
c51  =  Props (6) 

c61  =  Props (7) 

c52  =  Props (8) 

c62  =  Props (9) 

xnewdt  =  Props (10) 
inter  =  Props (11) 
relax  =  Props (12) 
toler  =  Props (13) 
c 
c 

c  Write  out  material  constants 
c 

write(NOUT,  '(///34h  Temperature  Independent  Constants  ,11)') 
write (NOUT, ' (15h  DNOT  ,  2x,  F10.3  )')  dO 

write (NOUT, ' (15h  Zl  ,  2x,  F10.3  )')  zl 

write (NOUT, ' (15h  C51  and  C61  ,  2(2x,  IPElO.l)')  C51,  C61 

write (NOUT, ' (15h  C52  and  C62  ,  2(2x,  IPElO.l)')  C52,  C62 

write (NOUT, ' (16h  Error  Tolerance,  lx,  IPElO.l)')  Toler 
write (NOUT, ' (16h  Relax.  Factor  ,  lx,  F10.3)')  relax 
write (NOUT, ' (16h  XNEWDT  ,  lx,  F10.3)')  XNEWDT 

c 

ipt  =  ntip  +  4 
c 

do  5  i  =  1,  NTSET 
T(i)  =  PROPS(ipt) 
ipt  =  ipt  +  1 

5  continue 
c 

do  6  i  =  1,  NTSET 
TE(i)  =  PROPS (ipt) 

ipt  =  ipt  +  1 

6  continue 
c 

do  7  i  =  1,  NTSET 
TNU(i)  =  PROPS (ipt) 
ipt  =  ipt  +  1 

7  continue 
c 

do  8  i  =  1,  NTSET 
TN(i)  =  PROPS(ipt) 

ipt  =  ipt  +  1 

8  continue 
c 

do  11  i  =  1,  NTSET 
TZ2(i)  =  PROPS (ipt) 
ipt  =  ipt  +  1 

11  continue 
c 

do  12  i  =  1,  NTSET 
TZ3(i)  =  PROPS (ipt) 
ipt  =  ipt  +  1 

12  continue 
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'  do  13  i  =  1-  NTSET 

TMl{i)  =  PROPS (ipt) 
ipt  =  ipt  +  1 

^3  continue 

'  do  16  i  =  1.  NTSET 

TM2(i)  =  PROPS (ipt) 
ipt  =  ipt  +  1 

1_6  continue 

"  do  18  i  =  1-  NTSET 

TRl(i)  =  PROPS (ipt) 
ipt  _  =  ipt  +  1 

18  continue 

""  do  19  i  =  1.  NTSET 

TR2(i)  =  PROPS (ipt) 
ipt  =  ipt  +  1 

19  continue 

^  do  20  i  =  1-  NTSET 

if  (inter  .gt.  0)  then 
TINTER(i)  =  inter 

TINTER(i)  =  PROPS  (rpt) 
ipt  =  ipt  +  1 

20  continue 

c  Write  Matrix  Elastic  Properties 

:  ilS  PoiSSns  R.tio  ,18h  spbcut  Factor 

”?a>  .TE<i)  .  TIMEEd, 

101  continue  „  >  , 

102  FORMAT(  3(F10.2,  8x) , 13x,  15) 

'l03  Boaaor-Parto.  Proportieo 

>  I8h  Temperature  . 12h  N 

>  12h  z3  ,12n  ri 

^  104  i  =  1.  NTSET 

^write(NOUT,105) 

104  continue  ^  o  u  t  > 

105  format ((  F10.2,  4x) , 6 (F12 . 3 ) ) 

\oo  Boaoer-Parto.  Propartija  Xl^ 

>  18h  Temperature  ,l/n 

^  do  107  i  =  1.  NTSET  _  ^ 

write (NOUT, 108)  T ( i) , TMl (I) ,TM2 (I) 

108  FORMAT((  F10.2,  4x) , 6 {F12 . 3 ) ) 

107  continue 
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return 

end 

c 

function  tprop(T,  PROPS,  Tget) 
c 

INCLUDE  'ABA_PARAM.INC' 
c 

parameter  (MAXP=20) 
real*8  T{MAXP),  PROPS (MAXP) 
c 

do  10  i  =  1,  (MAXP-1) 
if (tget  .eq.  t(i)  then 
tprop  =  props (i) 
return 

elseif  (tget  .eq.  t(i+l))  then 

tprop  =  props (i+1) 

return 

elseif  (tget  .gt.  t(i)  .and.  tget  .It.  t(i+l))  then, 
tprop  =  props (i) 

&  +  (PROPS (I+l) -PROPS (I) ) 

&  *  (Tget  -  T(I) ) 

&  /  (T(I+1)-T(I) ) 

return 
endif 

10  continue 

c 

write (6,*) 

write (6,*)  '  Temperature  Outside  of  Material  Data  Range  ' 

write (6,*) 
stop 
end 
c , 

subroutine  getprops (T1,DTEMP) 
c 

c 

c  Copyright  1995  The  University  of  Dayton 

c  All  rights  reserved, 

c 

c  This  material  may  be  reproduced  by  or  for  the 
c  U.S.  Govenment  pursuant  to  the  copyright  license 
c  under  the  clause  at  DFARS  252.227-7013  (Oct  1988) 
c 

c  Programmed  by : 
c 

c  Joseph  L.  Kroupa  (513)  255-9871 
c  University  of  Dayton  Research  Institute 
c  Dayton,  Ohio  45469-0128 
c 

c 

C  THIS  SOFTWARE  AND  ANY  ACCOMPANYING  DOCUMENTATION  IS  RELEASED  AS 

C  IS.  THE  U.S.  GOVERNMENT,  ITS  CONTRACTORS  AND  THEIR  SUBCONTRACTORS 

C  MAKE  NO  WARRANTY  OF  ANY  KIND,  EXPRESS  OR  IMPLIED,  CONCERNING  THIS 
C  SOFTWARE  AND  ANY  .ACCOMPANYING  DOCUMENTATION,  INCLUDING,  WITHOUT 
C  LIMITATION,  ANY  WARRANTIES  OF  MERCHANTABILITY  AND  FITNESS  FOR  A 
C  PARTICULAR  PURPOSE.  IN  NO  EVENT  WILL  THE  U.S.  GOVERNMENT,  ITS 
C  CONTRACTORS  AND  THEIR  SUBCONTRACTORS  BE  LIABLE  FOR  ANY  DAMAGES, 

C  INCLUDING  LOST  PROFITS,  LOST  SAVINGS  OR  OTHER  INCIDENTAL  OR 
C  CONSEQUENTIAL  DAMAGES  ARISING  OUT  OF  THE  USE,  OR  INABILITY  TO  USE 

C  THIS  SOFTWARE  OR  ANY  ACCOMPANYING  DOCUMENTATION,  EVEN  IF  INFORMED 
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C  IN  ADVANCE  OF  THE  POSSIBILITY  OF  SUCH  DAMAGES. 

C 

C  "This  software  is  being  used  at  the  user's  own  risk:  Neither  the 
C  Government  Agency  nor  its  contractors  assure  software's  accuracy 
C  or  its  appropriate  use . " 

C 


INCLUDE  'ABA_P ARAM. INC 
c 

PARAMETER  (MAXP=20) 

common  /matprop/  T(MAXP) ,  TE(MAXP),  TNU{MAXP) ,  TN(MAXP) , 

&  TZ2(MAXP).  TZ3(MAXP),  TMl(MAXP),  TM2 (MAXP) , 

&  TRl(MAXP),  TR2(MAXP),  TINTER ( MAXP ) , 

&  C51,  C61,  C52,  C62,  NTSET 

c 

COMMON  /matconst/  E,  G,  a3k,  anu,  dg,  d3k,  dnu,  dO, 

&  an,  ami,  am2,  zO,  zl,  z2,  z3 ,  al,  a2 ,  rl,  r2, 

&  dn,  dml,  dm2 , dzO , dz2 , dz3 , dal, da2 , drl , dr2 

c 

common  /control/  effhrd,  xnewdt,  relax,  toler,  inter,  Iplane 
c 

T2  =  T1  +  DTEMP 


c 

c 

c 


c 

c 

c 


C 

c 

c 


c 


c 

c 

c 


c 


c 


c 


c 


Extract  Elastic  Properties 

E  =  tprop ( T , TE , T1 ) 

E2  =  tprop (T,TE,T2) 
anu  =  tprop (T, TNU, Tl) 
anu2  =  tprop ( T , TNU , T2 ) 

G  =  E/(2.*(l.+  anu)-) 

note:  a3k  is  defined  as  3X  of  Bulk  Modulus 

a3k  =  E/(l.-2.*  anu) 

dg  =  E2 / ( 2 . * ( 1 . +  anu2 ) )  -  g 

d3k  =  E2/(l.-2.*  anu2)  -  a3k 

dnu  =  anu  -  anu2 

Special  Functional  curves  for  al  and  a2  defined  by  C5  and  C6 

al  =  c51  *  exp(-c61/ (Tl+273. ) ) 
dal  =  c51  *  exp(-c61/ (t2-t-273. )  )-al 

a2  =  c52  *  exp (-c62/ (Tl+273.)) 
da2  =  c52  *  exp (-c62/ (t2+273 . ) ) -a2 

Get  other  temperature  dependent  variables 


an 

= 

Tprop (T,TN, Tl) 

dn 

s: 

Tprop (T, TN, T2 ) 

-  an 

z2 

= 

Tprop (T,TZ2,T1) 

dz2 

=: 

Tprop(T,TZ2,T2) 

-  z2 

z3 

Tprop {T,TZ3,T1) 

dz3 

= 

Tprop (T,TZ3,T2) 

-  z3 

ami 

Tprop (T , TMl , Tl ) 

dml 

= 

Tprop { T , TMl , T2 ) 

-  ami 

am2 

= 

Tprop (T,TM2,T1) 
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dm2 

=  Tprop  {T,TM2,T2) 

-  am2 

rl 

=  Tprop {T,TR1,T1) 

drl 

=  Tprop (T,TR1,T2) 

-  rl 

r2 

=  Tprop (T,TR2,T1) 

dr2 

=  Tprop (T,TR2,T2) 

-  r2 

inter  =  1+IFIX (Tprop (T, TINTER, Tl) ) 
c 

return 

end 

c  end  of  user  subroutines 
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APPENDIX  B 


Listing  Of  Sample  Input  File 


♦HEADING 

Test  ir^Jut  file  for  urrat_db  -  BocSner  Partem,  imiaxial  tension 
♦mVEERCOT  MI13IMI2ATIOT,  SUPPRESS 
*N0DE,]SISErP=AU2I 


1 

0. 

0. 

0. 

2 

1. 

0. 

0. 

3 

1. 

1. 

0. 

4 

0. 

1. 

0. 

5 

0.5 

0. 

0. 

6 

1. 

0.5 

0. 

7 

0.5 

1. 

0. 

8 

0. 

0.5 

0. 

♦ELEMENT,  TYPE=CAX8 

,ELSE!r=ALLE 

1 

1  2 

3  4  5 

6 

♦SOLID  SECmCN,ELSET=ALLE,mTERIAL=B21S_3 

♦♦  Tiitetal21S  UDRI  Tref=900C 
♦♦  Special  Elastic  Case 


♦MATERIAL, 

NAME=B21S_0 

♦EXPANSICN, 

2ERD=900. 

9.770E-6 

23.0 

10.719E-6 

260.0 

11.500E-6 

482.0 

12.030E-6 

650.0 

12.348E-6 

760.0 

12.501E-6 

815.0 

12.271E-6 

900.0 

♦ELASTIC 

114000.0 

0.34 

23.0 

114000.0 

0.34 

260.0 

90000.0 

0.34 

482.0 

78000.0 

0.34 

650.0 

70000.0 

0.34 

760.0 

64000.0 

0.34 

815.0 

54700.0 

0.34 

900.0 

•kit 


♦♦  NEU'S  (1993)  TircEtal21S  CCNSTANTS 
** 


♦MATERIAL,  NAME=B21S_3 
♦♦♦EXPANSION,  ZERO=900.0 


** 

9.7787E-6 

23.0 

1.0713E-5 

260.0 

1.0915E-5 

315.0 

1.1093E-5 

365.0 

1.1267E-5 

415.0 

** 

1.1436E-5 

465.0 

★  ★ 

1.1492E-5 

482.0 

** 

1.1550E-5 

500.0 

1.1631E-5 

525.0 

1.1710E-5 

550.0 

★  ★ 

1.1788E-5 

575.0 

1.1865E-5 

600.0 

** 

1.2014E-5 

650.0 

★  * 

1.2323E-5 

760.0 
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**  1.2467E-5  815.0 

**  1.2689E-5  900.0 

*USER  MA.TEEIIAL,  CONSTANrS=176 

1  9  16  1.0E4 


1.37E4 

1.0 

23.0 

260.0 

525.0 

550.0 

112000.0 

108030.0 

95497.0 

93873.0 

0.340 

0.340 

0.340 

0.340 

4.80 

3.50 

1.280 

1.100 

1550.0 

1300.0 

1074.4 

1059.5 

100.0 

300.0 

1670.0 

2100.0 

0.0 

0.0 

0.0 

0.0 

0.350 

0.350 

6.822 

7.881 

3.0 

•  3.0 

3.0 

3.0 

3.0 

3.0 

3.0 

‘DEPVAR 

20 

*** 

3.0 

10 

0.75 

315.0 

365.0 

575.0 

600.0 

106130.0 

104090.0 

92172.0 

90395.0 

0.340 

0.340 

0.340 

0.340 

3.054 

2.649 

0.970 

0.820 

1250.4 

1205.4 

1044.6 

1029.8 

390.0 

500.0 

2600.0 

3700.0 

0.0 

0.0 

0.0 

0.0 

1.502 

2.549 

8.941 

10.00 

3.0 

3.0 

3.0 

3.0 

3.0 

3.0 

3.0 

3.0 

‘MATERIAL,  NAME=B21S_4 
★  ★ 

**  SANDER'S  (1995)  Tirnetal21S  CXONSTANTS 


‘EXPANSION, 

ZERD=900 

9.770E-6 

23.0 

10.719E-6 

260.0 

11.500E-6 

482.0 

12.030E-6 

650.0 

12.348E-6 

760.0 

12.501E-6 

815.0 

12.271E-6 

900.0 

‘USER  MATERIAL,  a)NSTANTS=76 

1.0E4 

1 

13 

6 

1.37E4 

2.0 

10 

0.75 

23.0 

260.0 

566.0 

650.0 

92700.0 

86612.0 

77216.0 

71964.0 

0.340 

0.340 

1.58 

1.35 

2950.0 

2650.0 

800.0 

200.0 

180.0 

180.0 

180.0 

180.0 

37.000 

50.000 

4.0 

4.0 

3.500 

3.500 

3.500 

3.500 

3.500 

3.500 

3.500 

3.500 

1600.0 

l.E-4 

415.0 

650.0 

101740.0 

86612.0 

0.340 

0.340 

2.243 

0.740 

1160.4 

1000.0 

660.0 

3800.0 

0.0 

0.0 

3.597 

10.00 

3.0 

3.0 

3.0 

.0 


3200.0 

l.E-4 

760.0 

0.340 

1.09 

190.0 

0.50 

4.0 

3.500 


5.8E5 

0.0 

465.0 

760.0 

99085.0 

77216.0 

0.340 

0.340 

1.838 

0.580 

1115.3 

600.0 

960.0 

4000.0 

0.0 

0.0 

4.644 

15.00 

3.0 

3.0 

3.0 

3.0 


1.37E4 
0.0 
482.0 
815.0 
98113.0 
71964.0 
0.340 
0.340 
1.700 
0.550 
1100.0 
300.0 
1100.0 
4100.0 
.  0.0 
0.0 
5.000 
30.00 
3.0 
3.0 
3.0 
3.0 


5.8E5 

0.0 

500.0 

900.0 

97045.0 

63122.0 

0.340 

0.340 

1.500 

0.550 

1089.3 

300.0 

1300.0 

4300.0 

0.0 

0.0 

5.763 

30.00 

3.0 

3.0 

3.0 

3.0 


5.5875E6 

0.00 

815.0 

0.340 

0.94 

150.0 

0.50 

4.0 

3.500 


1.37E4 

0.0 

112000.0 

0.340 

0.575 

180.0 

2.500 
4.0 

3.500 


5.5875E6 

0.0 

108030.0 

0.340 

0.55 

180.0 

7.000 

4.0 

3.500 


‘BOOSmARY 
1  1 

12 
2  2 

5  2 

4  1 

8  1 

***  Tsiperature  curve  for  case  1 

‘AMPLITUDE,  NAME=TaJRl 

0.00  25.0  48.0000  25.0 

★  ★★ 

‘“  Tarperatuore  curve  for  case  2  &  3 
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650.0 


*AMPLrTUDE,  NME=TCUR2 

0.00  650.0  4800.000 

**★ 

***  TsTpsTatuxe  cuive  for  case  4 
♦AMPLITUDE,  NAME=TCUR3 

0.00  25.00  24.0000 

72.000  25.00  96.0000 

***  Tenperature  curve  for  case  5 
♦AMPLITUDE,  NAME=TCUR4 

0.00  650.0  12.0000 
36.000  650.0  48.0000 

*** 

♦♦♦  Strain  curve  for  case  1,  2,  4, 

♦AMPLITUDE,  1©ME=MIUR1 

0.00  0.00  96.000 

♦♦♦  Strain  curve  for  case  3 
♦AMPLITUDE,  NME=MIUR2 

0.00  0.00  4800.000 

♦STEP,  MUNOTaNIC=ND,  INC=48 
This  is  a  tfcidirectional  Ttest  of  the  Directional  Bodner  Partem  Routines 

♦STATIC,  DIRECT 

1.0  48.0 

♦♦♦CONTROLS,  PARAMEIER=FIELD,  FIELD=GLOBAL 
♦*♦♦  0.0005 

♦BOUNDARY,  AMPLITUDE=MCUR1 

3  2  1.00 

4  2  1.00 

7  2  1.00 

♦TEMPERATURE,  AMPLITUDE==TCUR2 
ALLN  1.000 
♦EL  PRINT 
S 

SINV 

E 

♦NODE  PRINr,  ERBQUENCY=10 
U,RF 

♦FILE  FORMAT,  ASCII 

♦EL  FILE,  POSnTaN=AVERAGED  AT  NODES,  FREQUEWCY=1 

S,E,SINV,TEMP,  SDV 

♦NODE  FILE,  FREQUENCY=1 

COORD,  NT 

♦END  STEP 


25.00  36.00  482.0  60.000 

25.00 

650.0  18.00  760.0  30.000 

650.0 

&  5 

0.080 

0.040 


482.0 


760.0 
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APPENDIX  C 


Verification  Solution  of  Case  1 
(file  vcasel.math) 


Strain 

Stress  (MPa) 

Strain 

Stress  (MPa) 

0.0000000 

0.0000 

0.020832 

1145.6 

0.0008333 

93.330 

0.021666 

1145.9 

0.0016666 

186.66 

0.022499 

1146.2 

0.0024999 

279.99 

0.023332 

1146.3 

0.0033332 

373.32 

0.024166 

1146.4 

0.0041665 

466.65 

0.024999 

1146.5 

0.0049998 

559.98 

0.025832 

1146.6 

0.0058331 

653.31 

0.026666 

1146.6 

0.0066664 

746.64 

0.027499 

1146.7 

0.0074997 

839.97 

0.028332 

1146.7 

0.0083330 

933.30 

0.029166 

1146.7 

0.0091663 

1026.6 

0.029999 

1146.7 

0.0099996 

1083.9 

0.030832 

1146.7 

0.010833 

1098.6 

0.031665 

1146.7 

0.011666 

1110.4 

0.032499 

1146.7 

0.012500 

1119.6 

0.033332 

1146.7 

0.013333 

1126.7 

0.034165 

1146.7 

0.014166 

1132.0 

0.034999 

1146.7 

0.014999 

1136.0 

0.035832 

1146.7 

0.015833 

1138.9 

0.036665 

1146.7 

0.016666 

1141.1 

0.037499 

1146.7 

0.017499 

1142.7 

0.038332 

1146.7 

0.018333 

1143.8 

0.039165 

1146.7 

0.019166 

1144.6 

0.039998 

1146.7 

0.019999 

1145.2 
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APPENDIX  D 


Strain 

0.0000000 

0.0008333 

0.0016666 

0.0024999 

0.0033332 

0.0041665 

0.0049998 

0.0058331 

0.0066664 

0.0074997 

0.0083330 

0.0091663 

0.0099996 

0.010833 

0.011666 

0.012500 

0.013333 

0.014166 

0.014999 

0.015833 

0.016666 

0.017499 

0.018333 

0.019166 

0.019999 

0.020832 

0.021666 

0.022499 

0.023332 

0.024166 

0.024999 

0.025832 

0.026666 

0.027499 

0.028332 

0.029166 

0.029999 

0.030832 


Verification  Solution  Of  Case  2 


Stress  (MPa) 


Strain 


Stress  (MPa) 


0.0000 

0.031665 

383.51 

72.163 

0.032499 

383.51 

133.04 

0.033332 

383.51 

192.07 

0.034165 

383.51 

248.97 

0.034999 

383.51 

298.93 

0.035832 

383.51 

337.18 

0.036665 

383.51 

361.49 

0.037499 

383.51 

374.15 

0.038332 

383.51 

379.78 

0.039165 

383.51 

382.07 

0.039998 

383.51 

382.96 

383.30 

383.43 

383.48 

383.50 

383.51 
383.51 
383.51 
383.51 
383.51 
383.51 
383.51 
383.51 
383.51 
383.51 
383.51 
383.51 
383.51 
383.51 
383.51 
383.51 
383.51 
383.51 
383.51 
383.51 
383.51 
383.51 
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APPENDIX  E 


Verification  Solution  Of  Case  3 


Strain 

Stress  (MPa) 

0.0000000 

0.0000 

0.0008333 

72.106 

0.0016666 

121.44 

0.0024999 

138.00 

0.0033332 

140.14 

0.0041665 

140.36 

0.0049998 

140.38 

0.0058331 

140.38 

0.0066664 

140.38 

0.0074997 

140.38 

0.0083330 

140.38 

0.0091663 

140.38 

0.0099996 

140.38 

0.010833 

140.38 

0.011666 

140.38 

0.012500 

140.38 

0.013333 

140.38 

0.014166 

140.38 

0.014999 

140.38 

0.015833 

140.38 

0.016666 

140.38 

0.017499 

140.38 

0.018333 

140.38 

0.019166 

140.38 

0.019999 

140.38 

0.020832 

140.38 

0.021666 

140.38 

0.022499 

140.38 

0.023332 

140.38 

0.024166 

140.38 

0.024999 

140.38 

0.025832 

140.38 

0.026666 

140.38 

0.027499 

140.38 

0.028332 

140.38 

0.029166 

140.38 

0.029999 

140.38 

0.030832 

140.38 

Strain _ Stress  (MPa) 


0.031665 

140.38 

0.032499 

140.38 

0.034165 

140.38 

0.034999 

140.38 

0.035832 

140.38 

0.036665 

140.38 

0.037499 

140.38 

0.038332 

140.38 

0.039165 

140.38 

0.039998 

140.38 
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APPENDIX  F 


Verification  Solution  Of  Case  4 


Strain 

Stress 

(MPa) 

Temp  (°C) 

Strain 

Stress 

(MPa) 

Temp  (°C) 

0.00000000 

0.0000 

25.000 

0.030833 

787.05 

482.00 

0.00083333 

93.333 

25.000 

0.031667 

787.03 

482.00 

0.0016667 

186.67 

25.000 

0.032500 

787.03 

482.00 

0.0025000 

280.00 

25.000 

0.033333 

787.03 

482.00 

0.0033333 

373.33 

25.000 

0.034167 

787.03 

482.00 

0.0041667 

466.67 

25.000 

0.035000 

787.03 

482.00 

0.0050000 

560.00 

25.000 

0.035833 

787.03 

482.00 

0.0058333 

653.33 

25.000 

.  0.036667 

787.03 

482.00 

0.0066667 

746.67 

25.000 

0.037500 

787.03 

482.00 

0.0075000 

840.00 

25.000 

0.038333 

787.03 

482.00 

0.0083333 

933.33 

25.000 

0.039167 

787.03 

482.00 

0.0091667 

1026.7 

25.000 

0.040000 

787.03 

482.00 

0.0100000 

1083.9 

25.000 

0.040833 

787.03 

482.00 

0.010833 

1098.6 

25.000 

0.041667 

787.03 

482.00 

0.011667 

1110.4 

25.000 

0.042500 

787.03 

482.00 

0.012500 

1119.6 

25.000 

0.043333 

787.03 

482.00 

0.013333 

1126.7 

25.000 

0.044167 

787.03 

482.00 

0.014167 

1132.0 

25.000 

0.045000 

787.03 

482.00 

0.015000 

1136.0 

25.000 

0.045833 

787.03 

482.00 

0.015833 

1138.9 

25.000 

0.046667 

787.03 

482.00 

0.016667 

1141.1 

25.000 

0.047500 

787.03 

482.00 

0.017500 

1142.7 

25.000 

0.048333 

787.03 

482.00 

0.018333 

1143.8 

25.000 

0.049167 

787.03 

482.00 

0.019167 

1144.6 

25.000 

0.050000 

787.03 

482.00 

0.020000 

1145.2 

25.000 

0.050833 

821.57 

443.90 

0.020833 

1122.6 

63.100 

0.051667 

843.77 

405.80 

0.021667 

1097.2 

101.20 

0.052500 

879.05 

367.80 

0.022500 

1070.3 

139.30 

0.053333 

913.54 

329.70 

0.023333 

1041.9 

177.30 

0.054167 

945.70 

291.60 

0.024167 

1011.7 

215.40 

0.055000 

975.69 

253.50 

0.025000 

979.53 

253.50 

0.055833 

1008.1 

215.40 

0.025833 

949.64 

291.60 

0.056667 

1038.7 

177.30 

0.026667 

917.28 

329.70 

0.057500 

1067.4 

139.30 

0.027500 

882.53 

367.80 

0.058333 

1094.6 

101.20 

0.028333 

850.18 

405.80 

0.059167 

1120.3 

63.100 

0.029167 

825.88 

443.90 

0.060000 

1144.7 

25.000 

0.030000 

790.59 

482.00 

0.060833 

1146.7 

25.000 
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Strain 


Stress 

(MPa) 


0.061667 

1146.7 

0.062500 

1146.7 

0.063333 

1146.7 

0.064167 

1146.7 

0.065000 

1146.7 

0.065833 

1146.7 

0.066667 

1146.7 

0.067500 

1146.7 

0.068333 

1146.7 

0.069167 

1146.7 

0.070000 

1146.7 

0.070833 

1146.7 

0.071667 

1146.7 

0.072500 

1146.7 

0.073333 

1146.7 

0.074167 

1146.7 

0.075000 

1146.7 

0.075833 

1146.7 

0.076667 

1146.7 

0.077500 

1146.7 

0.078333 

1146.7 

0.079167 

1146.7 

0.080000 

1146.7 

Temp  (°C) 


25.000 

25.000 

25.000 

25.000 

25.000 

25.000 

25.000 

25.000 

25.000 

25.000 

25.000 

25.000 

25.000 

25.000 

25.000 

25.000 

25.000 

25.000 

25.000 

25.000 

25.000 

25.000 

25.000 
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APPENDIX  G 


Verification  Solution  Of  Case  5 


Strain 

Stress 

(MPa) 

Temp  (°C) 

Strain 

Stress 

(MPa) 

Temp  (°C) 

0.0000 

0.0000 

650.00 

0.015417 

136.98 

760.00 

0.00041667 

36.095 

650.00 

0.015833 

132.73 

760.00 

0.00083333 

72.190 

650.00 

0.016250 

130.64 

760.00 

0.0012500 

104.07 

650.00 

0.016667 

129.59 

760.00 

0.0016667 

133.08 

650.00 

0.017083 

129.05 

760.00 

0.0020833 

162.60 

650.00 

0.017500 

128.76 

760.00 

0.0025000 

192.14 

650.00 

0.017917 

128.62 

760.00 

0.0029167 

221.16 

650.00 

0.018333 

128.54 

760.00 

0.0033333 

249.09 

650.00 

0.018750 

128.50 

760.00 

0.0037500 

275.30 

650.00 

0.019167 

128.48 

760.00 

0.0041667 

299.15 

650.00 

0.019583 

128.47 

760.00 

0.0045833 

320.03 

650.00 

0.020000 

128.46 

760.00 

0.0050000 

337.52 

650.00 

0.020417 

128.46 

760.00 

0.0054167 

351.44 

650.00 

0.020833 

128.46 

760.00 

0.0058333 

361.94 

650.00 

0.021250 

128.46 

760.00 

0.0062500 

369.48 

650.00 

0.021667 

128.46 

760.00 

0.0066667 

374.65 

650.00 

0.022083 

128.46 

760.00 

0.0070833 

378.08 

650.00 

0.022500 

128.46 

760.00 

0.0075000 

380.30 

650.00 

0.022917 

128.46 

760.00 

0.0079167 

381.71 

650.00 

0.023333 

128.46 

760.00 

0.0083333 

382.60 

650.00 

0.023750 

128.46 

760.00 

0.0087500 

383.15 

650.00 

0.024167 

128.46 

760.00 

0.0091667 

383.49 

650.00 

0.024583 

128.46 

760.00 

0.0095833 

383.70 

650.00 

0.025000 

128.46 

760.00 

0.0100000 

383.83 

650.00 

0.025417 

134.59 

750.80 

0.010417 

374.19 

659.20 

0.025833 

144.80 

741.70 

0.010833 

357.33 

668.30 

0.026250 

157.42 

732.50 

0.011250 

336.61 

677.50 

0.026667 

171.74 

723.30 

0.011667 

314.21 

686.70 

0.027083 

187.43 

714.20 

0.012083 

291.34 

695.80 

0.027500 

204.29 

705.00 

0.012500 

268.60 

705.00 

0.027917 

222.20 

695.80 

0.012917 

246.28 

714.20 

0.028333 

241.07 

686.70 

0.013333 

224.60 

723.30 

0.028750 

260.85 

677.50 

0.013750 

203.67 

732.50 

0.029167 

281.50 

668.30 

0.014167 

183.61 

741.70 

0.029583 

302.96 

659.20 

0.014583 

164.50 

750.80 

0.030000 

325.21 

650.00 

0.015000 

146.43 

760.00 

0.030417 

341.88 

650.00 
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Strain 


Stress 


(MPa) 


0.030833 

354.79 

0.031250 

364.38 

0.031667 

371.18 

0.032083 

375.79 

0.032500 

378.83 

0.032917 

380.78 

0.033333 

382.01 

0.033750 

382.78 

0.034167 

383.26 

0.034583 

383.56 

0.035000 

383.74 

0.035417 

383.86 

0.035833 

383.93 

0.036250 

383.97 

0.036667 

383.99 

0.037083 

384.01 

0.037500 

384.02 

0.037917 

384.03 

0.038333 

384.03 

0.038750 

384.03 

0.039167 

384.04 

0.039583 

384.04 

0.040000 

384.04 

Temp  (°C) 


650.00 

650.00 

650.00 

650.00 

650.00 

650.00 

650.00 

650.00 

650.00 

650.00 

650.00 

650.00 

650.00 

650.00 

650.00 

650.00 

650.00 

650.00 

650.00 

650.00 

650.00 

650.00 

650.00 
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